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PREFACE 


This  program  is  intended  to  provide  a  completely  automatic 
process  for  reducing  “raw”  single-crystal  x-ray  diffraction 
intensity  data  to  the  form  of  structure  factors.  If  desired,  the 
data  reduction  may  be  extended  to  yield  normalized  structure  fac¬ 
tors.  The  program  is  applicable  primarily  to  data  obtained  on  film 
with  Weissenberg  and/or  precession  goniometers. 

All  programming  has  been  done  in  IBM  704  FORTRAN  II 
language.  Binary  program  decks  are  available  which  are  inde¬ 
pendent  of  any  special  operating  system. 

The  printed  output  has  been  contrived  to  be  easy  to  read  and 
to  give  the  program  user  easy  access  to  all  information  neces¬ 
sary  for  further  stages  of  the  structure  analysis.  Also  all  infor¬ 
mation  required  for  correction  of  the  initial  set  of  processed  data 
are  provided.  The  program  detects  some  errors  in  intensity  esti¬ 
mates  and  indexing,  and  indicates  these  both  in  a  preliminary 
error-detection  run  and  in  the  final  output. 

The  magnetic  tape  output  is  consistent  with  the  input  of  useful 
structure  determination  and  refinement  programs. 


PROBLEM  STATUS 

This  is  an  interim  report;  work  on  this  problem  is  continuing. 


AUTHORIZATION 

NRL  Problem  007 -03 
Project  RR  001-02-43-4805 


Manuscript  submitted  December  12,  1961. 
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AN  X-RAY  DIFFRACTION  DATA  REDUCTION 
PROGRAM  FOR  THE  IBM  704  AND  7000 


1.  CAPABILITIES 

Section  1  is  intended  to  give  a  general  impression  of  the  capabilities  at  the  program. 
Every  effort  has  been  made  to  make  the  program  as  general  as  possible,  and  it  is  felt 
that  all  contingencies  within  the  intended  scope  of  the  program  have  been  provided  for. 


1.1  Crystal  Settings 

Any  number  of  crystal  settings  can  be  accommodated  as  long  as  the  size  of  the  prob¬ 
lem  remains  within  the  limits  given  in  section  2.  The  goniometer  type,  x-ray  wavelength, 
and  zone  axis  indices  are  specified  for  each  setting.  General  settings  may  be  used;  e.g., 
a  setting  such  as  [  ll2]  may  be  specified  as  well  as  the  common  [  100],  [010],  [001] 
settings. 


1.2  Independently  Observed  Sets  of  Data 

Any  number  of  Independently  observed  sets  of  data  collected  from  a  given  crystal 
setting  may  be  accommodated  as  long  as  the  size  of  the  problem  remains  within  the  limits 
given  in  section  2. 


1.3  Intensity  Data 

The  program  accepts  ‘‘raw”  intensity  data  collected  by  multiple  film  Weissenberg  or 
precession  techniques.  That  is,  as  each  Intensity  estimate  is  read,  its  value  is  recorded 
on  a  data  sheet  which  will  be  given  to  a  key  punch  operator  without  further  manipulation 
by  the  crystallographer.  The  data  are  grouped  for  key  punching,  computer  input,  and 
processing  into  sets,  each  of  which  will  consist  of  data  read  from  a  pack  of  films  contain¬ 
ing  data  from  one  reciprocal  lattice  layer. 

The  data  from  the  different  films  in  the  multiple  film  packs  are  correlated  and  aver¬ 
aged.  The  correlation  is  done  Independently  for  each  pair  of  adjacent  films,  and  the  aver¬ 
aging  is  done  independently  for  each  set  of  data. 

Intensity  estimates  which  give  abnormal  ratios  of  intensities  (i.e.,  film  factors)  on 

adjacent  films  are  rejected  from  further  consideration  in  the  processing  unless  normal 
ratios,  involving  a  third  or  more  films,  also  are  found.  Rejection  or  acceptance  of  each 
film  factor  is  decided  by  comparison  of  each  film  factor  with  the  average  for  the  relevant 
pair  of  adjacent  films.  Acceptance  criteria  are  specified  by  the  program  user  for  each 
crystal  setting.  The  averaging  and  inspecting  process  is  iterated  until  no  further  rejec¬ 
tions  occur.  The  procedure  used  is  described  further  in  section  3.2.  Reflections  lying 
outside  the  range  of  observation  but  mistakenly  included  in  the  observed  set  also  are 
rejected.  Rejected  reflections  are  listed  in  a  separate  preliminary  output  especially 
designed  for  this  purpose,  and  also  are  clearly  labeled  in  the  final  data  output. 
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More  than  one  estimate  of  a  given  reflection  may  be  Included  in  any  set.  If  these  are 
to  be  averaged  later  (see  section  1.6)  they  must  be  given  identical  Miller  indices  (e.g., 
the  program  does  not  recognize  that  |F(hk{)|  =  |F(hk?)|,  see  section  2.1). 

Reflections  with  intensities  too  weak  to  be  observed  may  be  included.  Reflections 
absent  because  of  space  group  extinction  should  not  be  included. 

The  reflections  within  a  set  need  not  be  arranged  in  any  special  order  prior  to  load¬ 
ing  into  the  computer.  The  data  within  each  set  are  arranged  in  dictionary  order  of  h,  k,  -t 
in  the  first  phase  of  the  calculation.  However,  if  the  data  are  arranged  in  dictionary  order 
of  h,  k,  prior  to  loading  into  the  computer,  the  first  part  of  the  calculation  may  be  omitted. 

Within  any  independently  observed  set  or  "layer*  of  data,  two  different  intensity  scales 
may  be  used  for  estimating  Intensities.  For  example,  one  scale  may  cover  the  range  of 
intense  reflections  while  the  other  scale  covers  the  range  of  weak  reflections.  The  corre¬ 
lation  factor  relating  the  two  scales  must  be  known. 

Any  independently  observed  set  of  intensity  data  may  be  on  a  arbitrary  scale  which 
has  not  been  correlated  with  any  other  set  of  data. 

The  data  may  be  loaded  into  the  computer  from  either  punched  cards  or  magnetic  tape. 


1.4  Lp  and  Spot  Shape  Corrections 

The  data  are  corrected  for  the  Lorentz  and  polarization  effects.  For  Weissenberg 
data  correction  may  or  may  not  be  made  for  compaction  or  expansion  of  spot  shapes  on 
upper  layer  photographs.  Where  both  expanded  and  compacted  spots  are  estimated  from 
the  same  films,  the  compacted  data  must  be  treated  as  a  set  separate  from  the  expanded 
data.  The  method  of  Phillips  (1)  is  used  for  this  correction.  Where  possible,  the  expanded 
spots  should  be  used  because  the  correction  factor  for  compacted  spots  passes  through 
zero  at  low  (sin  e)/K.  When  the  spot  shape  correction  factor  is  equal  to  or  less  'um  10*’, 
the  Intensity  for  that  reflection  is  set  equal  to  zero  and  the  reflection  is  listed  witn  the 
-ejected  data. 

The  equations  of  Waser  (2)  are  used  for  the  precession  Lp  Correctims. 


1.5  Absorption  Corrections 

The  program  does  not  calculate  absorption  corrections,  but  provision  is  made  in  the 
input  for  inclusion  of  previously  determined  absorption  corrections.  The  corrections  are 
applied  automatically  in  the  computer. 


1.6  Correlation  of  Independently  Observed  Sets  of  Data 

if  more  than  one  crystal  setting  has  been  used  in  collecting  the  data,  the  independently 
observed  sets  of  data  are  correlated  and  put  on  a  common  arbitrary  scale  by  the  least 
squares  procedure  of  Rollett  and  Sparks  (3).  The  weighting  of  observations  is  done  as 
follows:,  The  relevant  observations  consist  of  pairs  of  nonzero  Intensities  having  identical 
sets  of  'Miller  Indices,  each  intensity  of  a  pair  having  been  obtained  from  a  different  crystal 
setting.  i2ach  uitensiiy  of  a  pair  is  multiplied  by  the  factor  (n,  +  +  1,,),  where  the 

subscripts  a  and  b  identify  the  individuals  of  the  pair,  n  is  the  number  of  films  from  iHiich 
the  intensity  has  been  estimated,  and  i  is  the  value  of  the  intensity  estimate  corrected 
for  Lp  and  spot  shape.  The  reflection  data  must  be  arranged  in  dictionary  order  of  h,  k, 
at  this  point  in  the  calculation.  The  eigenvalue -eigenvector  solution  is  obtained  by  a 
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modified  version  of  SHARE  program  AN  F  202  (distribution  number  664)  and  SHARE  pro 
gram  RW  VCTR  (distribution  number  635).  These  programs  are  FORTRAN  coded. 


1.7  Final  Averaging  of  Data 

After  having  been  scaled  to  a  common  basis,  all  equivalent  data  (i.e.,  all  data  with 
identical  Miller  index  triples)  are  averaged  and  the  data  are  arranged  in  a  single  listing 
in  dictionary  order  of  h,  k,  -t.  In  computing  average  intensities,  the  individual  intensities 
are  weighted  according  to  the  number  of  films  from  which  each  was  estimated.  The 
data  also  are  maintained  in  unaveraged  form  in  the  sets  or  '‘layers”  as  originally 
observed. 


1.8  Placing  the  Data  on  the  Absolute  Scale 

The  data  are  placed  on  the  absolute  scale  by  means  of  a  modified  Wilson  function 
suggested  by  Schoomer  and  Hauptman  (4).  The  function  used  is 

loR  Kj  =  A  +  BSj’* 

where  A,  B,andX  are  constants  determined  as  described  below,  is  the  value  of  (sin  d)^/k^ 
at  the  middle  of  the  ith  increment  in  (sin  0)  Vk}  within  which  values  of  are  calculated, 
and 


where  c  is  the  integer  factor  used  to  convert  (the  normalized  structure  factor  (5))  to  8^ 
(the  quasi -normalized  structure  factor  (6)),  i  is  the  value  on  an  arbitrary  scale,  and 
cr,  =  if?.  The  summations  are  taken  over  all  reflections  observed  in  the  ith  increment 

of  (sin  I") 

The  constants  A  and  n  are  found  by  calculating  least  squares  fits  lor  values  of  X 
ranging  from  1  through  10  in  increments  of  0.5.  The  “best  parameters”  are  those  yielding 
the  smallest  sum  of  squares  of  residuals.  In  the  least  squares  calculations  each  obser¬ 
vation  is  weighted  according  to  the  number  of  reflections  contributing  to  it. 

A  fixed  increment  in  (sin  0)'^/k^  is  used  in  the  calculations.  This  increment  is  speci¬ 
fied  by  the  program  user.  The  range  of  the  function  is  to  the  edge  of  the  sphere  of  reflec¬ 
tion  implied  by  the  smallest  wavelength  used.  Reflection  data  need  not  be  observed  in  all 

increments  in  this  range. 

The  atomic  scattering  factors  f  are  calculated  using  the  equation 

f  A  exp  [-n  (sin +  B  exp  [“b  (s  i  n  +C 

proposed  by  Vand,  Eiland,  and  Pepinsky(7).  Values  of  the  constants  A,  .n,  b,  b,  and  C  have 
been  calculated  and  tabulated  for  a  large  number  of  atoms  and  ions  by  Forsyth  and  Wells 
(8).  The  approximations  are  good  through  the  molybdenum  range  of  (sin  0)/k. 


1.9  Preliminary  Output 

The  preliminary  output  occurs  during  the  phase  of  the  calculation  consisting  of  filni- 
to-film  averaging  of  intensity  data  and  correction  for  Lp  and  spot  shape  effects.  This 
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output  may  be  observed  in  a  pilot  run  for  the  purpose  of  error  detection.  In  this  case, 
the  data  do  not  need  to  be  arranged  in  dictionary  order  of  h,  k,  I,  and  the  first  part  of  the 
calculation,  DICNRY,  may  be  omitted.  The  data  processing  may  be  stopped  by  sense 
switch  control  after  completion  of  the  preliminary  output.  This  output  consists  of: 

A.  The  film  factors  for  each  independent  set  of  data  (see  section  1.3).  These  are 
printed  immediately  after  the  Lp  corrections  have  been  applied  and  before  the  next  set  of 
data  are  read  into  the  computer.  The  output  for  each  layer  or  set  of  data  consists  of: 

1.  the  zone  axis  indices,  u,  V,  W 

2.  the  layer  number,  LN 

3.  the  number  of  reflections  in  the  set,  NR 

4.  the  five  film  factors. 

B.  Rejected  reflections  are  listed,  along  with  sufficient  information  to  identify  their 
original  sets,  after  all  of  the  data  have  been  processed  as  far  as  application  of  the  Lp 
corrections.  The  output  for  each  layer  or  set  of  data  consists  of: 

1.  U,  V,  w 

2.  LN. 

The  output  for  each  rejected  reflection  is: 

1.  h,  k,  -i 

2.  (sin  6) 

3.  the  square  root  of  the  Lp  correction  factor  (including  the  spot  shape  correction  factor) 

4.  IfI* 

5.  A  column  headed  with  R  in  which  rejections  are  indicated  as  follows: 

a.  A  1  indicates  that  the  rejection  is  a  result  of  an  abnormal  film  factor. 

b.  A  10  indicates,  for  Weissenberg  data,  that  the  reflection  lies  outside  the  sphere 
of  reflection  implied  by  the  x-ray  wavelength  specified  for  the  relevant  crystal  setting. 

For  precession  data,  a  10  indicates  that  the  reflection  lies  beyond  the  range  of  accessi¬ 
bility  of  the  goniometer  for  the  particular  settings  specified. 

c.  A  20  indicates  that  the  spot  shape  correction  factor  is  equal  to  or  less  than 
10"^  This  type  of  rejection  can  occur  only  for  compacted  reflections  observed  on  upper 
layer  Weissenberg  photographs. 

The  rejection  numbers  are  additive.  Thus  an  11  would  indicate  that  the  reflection 
datum  suffers  from  both  an  anomalous  film  factor  and  incorrect  indexing.  A  20  could 
possibly  indicate  that  the  relfection  has  been  observed  twice  with  incorrect  indexing, 
etc.  In  any  case,  the  unaveraged  data  listed  in  the  sets,  as  originally  observed  (see 
section  I.IO.C)  should  be  consulted  to  determine  the  exact  nature  and  location  of  the 
difficulty. 

The  preliminary  output  may  be  observed  on-line  if  so  desired.  In  any  case,  it  is 
written  on  the  BCD  output  tape  for  off-line  printing. 


NAVAL  RISiAftCH  LABORATORY 


5 


1.10  Output  of  Structure  Factor  Data 

The  BCD  output  of  structure  factor  data  for  off-line  printing  nnay  be  divided  Into  three 
classifications,  A,  B,  C,  as  follows: 

A.  The  data  are  listed  in  dictionary  order  of  h,  k,  -f.  All  |f|^  data  with  equivalent 

h,  k,  i  are  averaged  in  this  output  (see  section  1.7).  Each  independent  reflection  is  listed 
only  once.  The  output  for  each  reflection  is: 

1.  h,  k,  -t 

2.  (sin  6f)/\ 

3.  !fP 

4.  |f1 

5.  t,  i.e.,  the  integer  factor  required  to  convert  E*  to  £^. 

6.  The  total  number  of  filmsfrom  which  nonzero  estimates  of  the  intensity  were  made. 
This  column  is  headed  w,  and  these  quantities  will  be  hereafter  designated  by  w. 

Following  the  last  reflection  in  this  listing  is  the  total  number,  modulo  32728,  of 
independent  data  printed  (including  rejected  data). 

B.  The  modified  Wilson  equation  parameters  and  sufficient  data  to  hand-plot  the 
curve”  -  specifically: 

1.  A,  B,  X  from  the  equation 


Log  K  =  A  +  B  ((sin 

2.  The  “absolute  base”  scale  factor  determined  from  A,  i.e.,  the  factor  by  which 
the  correlated  |fP  data  are  to  be  multiplied  to  put  them  on  the  absolute  scale. 

3.  Midincrement  values  of  (sin  oy/K? . 

4.  The  observed  values  of  tog  k  (Naperian  base). 

5.  The  numbers  of  reflections  observed  within  each  increment  of  (sin  0)*A*. 

C.  Finally,  the  data  are  listed  in  the  sets  or  “layers”  of  Independently  observed  data 
as  they  were  fed  into  the  computer.  These  data  have  not  been  averaged  with  other  equiva¬ 
lent  data.  The  data  within  each  set  are  listed  in  dictionary  order  of  h,  k,  At  the  head 
of  each  set  are  the  zone  indices  [uvw],  the  layer  number,  the  printed  words  PRECESSION 
or  WEISSENBERG  to  identify  the  goniometer  type,  the  number  of  data  in  the  set,  and  a 
scale  factor  which  is  the  product  of  the  layer  least  squares  correlation  factor  and  the 
“absolute  base”  scale  factor  determined  from  the  modified  Wilson  function.  The  output  for 
each  reflection  is: 

1.  Ii,  k,  't 

2.  (sin  0)/^ 

3.  The  square  root  of  the  Lp  factor  (including  the  spot  shape  correction  factor  if 
this  has  been  applied). 

4.  iFl^ 

5.  1f1 

6.  f 
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7.  w  and  R  as  explained  previously  except  that  these  quantities  now  are  applicable 
to  the  individual  sets  of  data  only  (see  sections  1.10. A  and  1.9.B). 

In  addition  to  the  printed  output,  two  BCD  tapes  are  written  in  a  format  consistent 
with  the  least  squares  refinement  program  of  Busing  and  Levy  (9).  These  tapes  contain: 

Tape  5  -  The  av^aged  data  in  the  form  of  h,  v.,  I,  IfI^,  a  least  squares  weighting  factor 
which  is  equal  to  |f|Aw,  where  w  is  defined  above  (section  I.IO.A),  and(sin  e)/k  (incolunms 
55-63  to  be  stored  in  T 13  (9)).  The  data  are  listed  in  the  same  order  as  in  the  printed  output. 

Tape  7  -  The  same  Information  for  the  unaveraged  independent  sets  of  data.  The 
data  are  listed  in  the  same  order  as  in  the  printed  output. 


1.11  Normalization  of  the  Structure  Factors 

If  so  desired,  the  data  reduction  may  be  extended  to  yield  functions  of  |e|.  The  con¬ 
version  from  |f|  values  to  |E|  values  is  done  as  follows:  The  |f|^  data  are  corrected  for 
the  temperature  effect  using  the  analytical  expression  for  the  K  curve  as  ^plied  by  the 
equations  in  section  1.8.  The  resulting  |f|2  values  are  then  divided  by  e  1  f?.  The  values 
of  f  are  calculated  as  described  in  section  1.8. 

1.12  Output  of  Normalized  Structure  Factor  Data 

The  printed  output  is  quite  similar  to  the  printed  output  of  f  data  already  described. 

It  is  divided  into  three  parts,  A,  B,  and  C,  as  follows: 

A.  The  averaged  data  are  listed  in  dictionary  order  of  h,  k,  Bach  independent 
reflection  is  listed  only  once.  The  output  for  each  reflection  is: 

1.  h,  k,  i 

2.  (sin  0)/X 

3.  I  El*  -  1 

4.  |e1 

5.  e ,  w,  and  R  as  defined  previously  (see  section  I.IO.A). 

B.  The  “K  curve”  parameters  and  average  values  for  several  functions  of  |e|.* 
Specifically: 

1.  A,  B,  X  from  the  equation  toKK  =  A  +  b[ (sin  d)^/K^]’' 

2.  The  “absolute  base”  scale  factor  determined  from  A. 

3.  <|E|> 

4.  <I|EI*-11> 

5.  <|EI*-1> 

C.  Finally,  the  data  are  listed  in  the  sets  or  “layers”  of  independently  observed  data 
as  they  were  fed  into  the  computer.  The  data  within  each  set  are  listed  in  dictionary  order 
of  h,  k,  -i.  These  data  have  not  been  averaged  with  other  equivalent  data.  At  the  head  of 


*The  program  does  not  take  cognizance  of  the  fact  that  for  structures  giving  rise  to  both 
real  and  complex  structure  factors  the  averages  for  the  corresponding  sets  of  intensities 
are  different.  Usually,  the  proportion  of  real  structure  factors  is  small  and  this  factor 
may  be  ignored;  however,  occasionally  it  will  have  to  be  taken  into  account. 
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each  set  are  the  zone  indices  fuvw],  the  layer  number,  the  printed  words  PRECESSION 
or  WEISSENBERG  to  identify  the  goniometer  type,  and  the  number  of  reflections  in  the 
set.  The  output  for  each  reflection  is: 

1.  h,  k, 

2.  (s  ill  0)/k 

3.  A  factor  K  by  which  the  |Fp  values,  already  scaled  on  the  absolute  base,  are  to 
be  multiplied  to  convert  them  to  |e1^  values. 

4.  lEl^ -  1 

5.  If. I 

6.  e,  w,  R  as  defined  previously  (see  section  1.10. A). 

In  addition  to  the  printed  output,  two  BCD  tapes  are  written  in  the  format  of  the 
Fourier  program  of  Sly  and  Shoemaker  (10).  These  tapes  contain: 

Tape  5  -  The  averaged  data,  listed  in  the  same  order  as  in  the  printed  output,  in  the 
form  of  h,  k,  i  and  |Ep  -  1. 

Tape  7  -  The  same  information  written  on  tape  for  the  unaveraged  independent  sets 
of  data  in  the  same  order  as  in  the  printed  output. 


1.13  Ordering  of  |E|  and  |eP-  1  Data 

If  so  desired,  the  data  reduction  may  be  extended  to  an  arranging  and  output  of  the 
data  in  decreasing  order  of  magnitude  of  |eP  -  1.  These  data  are  then  segregMed 
into  groups  specified  according  to  evenness  («)  or  oddness  (u)  of  h,  k,  and  C.  Omy  data 
for  which  |e|^  -  1  5-0  are  listed. 

The  output  consists  of  the  following  groupings: 


1. 

All  data 

4. 

Ki'R 

data 

7. 

URll 

data 

2. 

Kki;  data 

5. 

KRU 

data 

8. 

RUll 

data 

3. 

uKK  data 

6. 

IIIIR 

data 

9. 

null 

data. 

In  all  groups,  the  data  are  listed  in  decreasing  order  of  magnitude  of  |e|^  -  1. 


1.14  Time  of  Computation 

No  formulas  for  estimating  computing  time  are  available.  However,  the  following 
examples  should  suffice  to  give  a  general  impression.  The  times  given  include  time  used 
in  card  program  read-in  and  tape  changes  between  the  different  parts  of  the  calculations. 
The  data  were  processed  completely,  from  DICNRY  through  E  ORDER. 


Crystal 

Settings 

Sets 

Input 

Reflections 

Time 

Computer 

Cyclohexaglycyl 

1 

1 

407 

20  min 

704 

6-Hydroxycrinamine 

2 

11 

4195 

90  min 

704 

6-Hydroxycrinamine 

2 

11 

4195 

30  min 

7090 
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2.  LIMITATIONS 

2.1  Indexing  of  Reflection  Data 

Care  must  be  taken  to  maintain  consistency  in  the  assignment  of  Miller  indices.  The 
computer  receives  no  crystallographic  symmetry  information.  Therefore,  crystallograph- 
Ically  equivalent  reflections  must  be  assigned  identical  index  triples  to  assure  proper 
operation  of  the  layer  correlation  and  final  intensity  averaging  calculations.  Furthermore, 
the  Miller  index  assignments  must  be  consistent  with  the  zone  axis  indices  in  a  certain 
manner.  The  criterion  of  consistency  is  that 

bU  +  kV  +  4W  5-0 

where  h,  k,  -I  are  the  Miller  indices,  and  u,  v,  w  are  the  zone  axis  indices  for  the  relevant 
crystal  setting. 


2.2  Data  Reduction  to  Structure  Factor  Form 

1.  The  processing  is  limited  to  data  collected  by  goniometers  having  the  Weissenberg 
or  precession  geometries.  Though  primarily  designed  for  processing  film  data,  it  will 
accommodate  counter  data  obtained  with  Weissenberg  or  precession  goniometers  (e.g., 
the  Geiger  counter  goniometer  described  by  Evans  (11)). 

2.  The  program  does  not  provide  for  calculation  of  absorption  corrections.  However, 
previously  determined  absorption  corrections  may  be  included  in  the  input  data  for  each 
reflection,  and  these  corrections  are  automatically  applied. 

3.  ijtc  number  of  different  atomic  species  in  the  crystal  must  hot  exceed  20. 

4.  No  more  than  ten  different  x-radiation  wavelengths  can  be  used  in  obtaining  the 
data. 

5.  Number  of  data:  Each  independent  set  or  “layer”  of  data  may  contain  as  many 

as  iOOO  reflections.  Sets  containing  more  than  1000  reflections  may  be  divided  into  smaller 
sets  of  1000  or  less  reflections.  There  may  be  as  many  as  50  sets  of  data.  In  total,  the 
program  will  accommodate  as  many  as  50,000  reflections  involving  as  many  as  300,000 
intensity  estimates. 

6.  The  number  of  films  from  which  any  set  of  intensities  has  been  estimated  cannot 
exceed  six. 

7.  If  there  is  more  than  one  crystal  setting,  each  set  of  data  must  contain  at  least 
one  nonzero  reflection  which  also  is  observed  with  nonzero  value  from  another  crystal 
setting.  The  number  of  input  reflections  (as  distinct  from  unique  reflections)  in  common 
between  two  reciprocal  lattice  nets  from  different  settings  must  not  exceed  100  for  each 
net.* 

8.  When  one  or  several  sets  of  data  obviously  are  off-scale  from  the  bulk  of  the  data 
by  as  much  as  two  or  more  orders  of  magnitude  (as  could  happen  when  very  short  expo¬ 
sures  are  correlated  with  very  long  exposures)  a  preliminary  scaling,  sufficient  to  scale 
the  data  to  rouglily  the  same  order  of  magnitude,  is  advisable.  This  may  be  done  very 
conveniently  by  including  the  scale  factor  in  the  absorption  correction  factor. 


’•'When  a  crystal  setting  is  repeated  in  the  input,  as  would  happen  when  one  or  several 
layers  are  observed  with  both  Weissenberg  and  precession  goniometers,  cross - 
correlation  of  these  particular  layers  is  not  attempted  by  the  program.  However,  all 
other  correlation  data  involving  these  layers  is  used  in  the  usual  manner. 
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9.  The  total  number  of  input  reflections  having  the  same  Miller  index  h  in  common 
must  not  exceed  1000.  The  total  number  of  input  reflections  having  both  h  and  l:  in  com¬ 
mon  must  not  exceed  100. 

10.  The  exponent  X  in  the  modified  Wilson  equation 

log  K  =  A  B  [(sin  * 

is  restricted  to  integral  and  half  integral  values  between  1  and  10  inclusive. 

11.  In  evaluating  the  parameters  A,  B,  and  X  in  the  modified  Wilson  equation,  the 

number  of  intervals  in  (sin  within  each  of  which  K  is  evaluated  from  the  experi¬ 

mental  data,  must  not  exceed  100. 

12.  Data  obtained  with  radiation  of  wavelength  shorter  than  that  of  MoKa  may  be 
processed,  but  those  parts  of  the  data  processing  requiring  evaluation  of  atomic  scattering 
factors  will  have  dubious  reliability  because  of  the  limitations  of  the  Vand,  Eiland, 
Pepinsky  approximation  equation  for  calculation  of  these  quantities. 


2.3  Data  Reduction  to  the  Form  of  Normalized  Structure  Factors 
No  limitations  in  addition  to  those  above  apply. 


2.4  Data  Reduction  Extended  to  Include  Ordering  of  |e|^  -  1 
The  following  limitations  apply  in  addition  to  those  above: 

1.  A  maximum  of  5000  unique  reflections  (including  reflections  for  which  |e|*  - 
1  <  0)  can  be  accommodated  in  one  pass.  Of  these  5000  reflections  no  more  than  2400 
can  have  |e|*  -  1  s-  0.  Any  number  of  passes  may  be  made.  Each  pass  orders  2400,  or 
less,  reflections  independently  of  the  reflections  ordered  in  other  passes.  If  there  are 
more  than  5000  reflections  on  the  input  tape,  the  program  automatically  orders  the  data 
in  groups  of  5000  (or  less  for  the  last  group)  until  all  of  the  data  are  ordered. 

2.  No  more  than  300  reflections  can  occur  in  any  one  of  the  eigdrt  groups  segregated 
according  to  evenness  or  oddness  of  Miller  Indices  for  any  one  pass. 


3.  INPUT  DATA  AND  FORMATS 

3.1  Three  Categories  of  Input  Data 

The  input  data  are  categorized  into  three  main  parts,  each  succeeding  part  being  a 
subdivision  of  the  preceding  one,  as  follows: 


I.  General  Information 

General  crystallographic  and  chemical  information. 

II.  Crystal  Setting  Information 

1.  Information  which  defines  an  individual  crystal  setting. 

2.  Information  which  is  pertinent,  in  general,  only  to  all  data  collected  from  this  setting. 
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III.  Reciprocal  Net  Information 

1.  Information  which  defines  an  individual  set  of  reflection  data.  This  set  of  data 
consists  of  all  of  or  part  of  the  observable  portion  of  a  reciprecal  lattice  net  obtained  at 
the  crystal  setting  specified  in  II. 

2.  Information  which  is  pertinent  in  general,  only  to  this  particular  set  of  data. 

3.  The  reflection  data. 

The  program  user  should  understand  clearly  that  the  data  in  category  I  are  specified 
only  once  for  a  given  crystal.  Category  II  data  are  specified  for  each  crystal  setting,  and 
category  III  data  are  specified  for  each  individual  set  or  “layer”  of  data  obtained  at  each 
setting. 

In  the  exposition  of  the  three  categories  given  below,  each  type  of  card  to  be  punched 
is  treated  in  detail  and  is  identified  by  number.  A  diagrammatic  presentation  of  card 
composition  with  FORTRAN  formats  is  Included;  the  sequence  of  cards  in  the  data  deck 
is  given  last. 


3.2  Category  I  -  General  Information 
Card  1: 

Card  1  is  for  input  of  Hollerith  information,  e.g.,  the  name  of  the  crystal  and  the 
name  of  the  crystallographcr.  A  maximum  of  72  characters  are  allowed. 


Card  2: 


columns 

1  -  11 

12  -  22 

23  -  33 

34  -  44 

45  -  55 

56  -  66 

data 

At 

A2 

A3 

al 

a2 

a3 

format 

F11.7 

F11.7 

Fll.7 

Fll.7 

Fll.7 

Fll.7 

Card  2  contains  the  “real-space”  unit  cell  parameters  in  the  order  a,  b,  c,  a,  y. 

The  parameters  a,  b,  and  c  are  in  units  of  angstroms;  a,  fs,  and  y  are  expressed  in  degrees 
and  are  always  given,  even  when  they  are  fixed  by  space  group  symmetry. 


columns 

1  -  12 

13  -  24 

data 

R1 

R2 

format 

E12.5 

E12.5 

Card  3  contains  the  x-ray  beam  collimator  parameters  R1  and  R2  as  defined  by 
Phillips  (1)  for  use  in  correcting  intensity  data  for  the  expansion  or  compaction  of  spot 
size  on  upper  layer  equi-inclination  Weissenberg  photographs.  For  the  standard 
Weissenberg  goniometer  these  parameters  have  the  values  R1  =  28.648  and  R2  =  78.000. 
Where  this  correction  is  not  to  be  applied  or  is  irrelevant,  a  blank  card  may  be  used. 
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Card  4: 


columns 

1  -  4 

5-9 

data 

NS 

S 

format 

14 

F5.4 

Card  4  contains: 

1.  The  number  of  different  atomic  species  present  in  the  crystal,  NS. 

2.  The  interval  in  sin^  s,  to  be  used  in  calculating  the  log  K  curve  for  deter¬ 
mination  of  the  absolute  base  acale  factor.  The  unit  is  angstrom'^. 


Card  5: 


columns 

1  -  8 

9-16 

17  -  24 

25  -  32 

33-40 

41  -  44 

data 

A 

a 

B 

b 

C 

NA 

format 

F8.4 

F8.4 

F8.4 

F8.4 

F8.4 

F4.0 

Card  5  is  applicable  to  a  single  atomic  species.  There  is  a  separate  card  for  each 
species.  The  card  contains: 

1.  The  parameters  A,  a,  b,  b,  and  c  as  defined  by  Forsyth  and  Wells  (8)for  calcula¬ 
tion  of  the  atomic  scattering  factor.  Values  of  these  parameters  are  tabulated  in  the 
paper  by  Forsyth  and  Wells  (see  section  1.8). 

2.  The  number  of  atoms  of  this  atomic  species  in  the  unit  cell,  NA. 


3.3  Category  II  -  Crystal  Setting  Information 
Card  6: 


columns 

1  -  3 

4-6 

7-9 

10  -  12 

13  -  15 

16  -  26 

27  -  37 

38-48 

data 

U 

V 

W 

ST 

GT 

RF 

L 

RS 

format 

F3.0 

F3.0 

F3.0 

F3.0 

F3.0 

F11.8 

F11.8 

F11.8 

Card  6  contains: 

1.  The  zone  axis  indices  u,  V,  w.  As  an  example,  consider  a  crystal  mounted  on 
its  a  axis.  The  zone  axis  indices  for  Weissenberg  data  would  be  1,  0,  0.  For  the  same 
crystal  aligned  on  a  precession  goniometer  so  as  to  obtain  the  ho^,  bi^,  h2l,  etc.  nets, 
the  zone  axis  indices  would  be  0,  1,  0  (see  section  2.1). 

2.  A  test  flag,  ST,  to  answer  the  question:  Is  this  the  last  setting  to  be  read  into 
the  computer? 

Yes:  ST  =  1 
No:  leave  blank. 

3.  A  test  flag,  GT,  to  answer  the  question:  Is  this  Weissenberg  or  precession  data? 

precession:  CT  1 

Weissenberg:  leave  blank. 
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4.  The  rejection  factor,  RF.  In  determing  the  average  film  factor  (i.e.,  ratio  of 
intensity  estimates)  between  a^acent  films  in  a  pack,  the  ratios  of  all  individual  nonzero 
intensities  appearing  on  both  films  are  calculated  and  averaged.  The  product  of  this  aver¬ 
age  and  the  “rejection  factor”  is  the  upper  limit  of  acceptance  of  the  ratios.  The  quotient 
of  the  average  and  the  “rejection  factor”  is  the  lower  limit  of  acceptance.  Both  intensity 
readings  of  a  pair  which  yield  a  film  factor  falling  outside  the  range  of  acceptance  are 
rejected  unless  one  of  them  also  yields  an  .acceptable  film  factor  involving  a  third  film. 

See  section  1.3  for  additional  details. 

5.  Wavelength  of  x-radiation  used,  L,  expressed  in  Angstrom  units. 

6.  A  value  of  (ain  d)/K,  RS,  is  specified  at  and  above  which  the  aj  and  aj  spots 
are  clearly  resolved.  Intensities  of  reflections  with  (ain  e)/K  equal  to  or  greater  than 
this  value  are  multiplied  by  3/2.  The  unit  is  angstrom* 


Card  7: 


colunms 

1  -  3 

4-6 

7-9 

10  -  20 

21  -  31 

data 

e 

f 

a 

9 

format 

F3.0 

F3.0 

F3.0 

F11.7 

F11.7 

Card  7  is  included  if  and  only  if  the  data  are  obtained  by  the  precession  method.  The 
card  contains: 

1.  The  net  orientation  Indices  e,  f,  g,  which  are  the  Miller  indices  of  the  lowest 
order  reciprocal  lattice  point  (even  though  its  reflection  may  not  appear  on  the  film  for 
any  reason)  lying  on  the  upper  half  of  the  most  nearly  vertical  central  lattice  line  appear¬ 
ing  in  the  zero  layer  reciprocal  net  (see  section  2.1).  This  line  is  designated  as  the  ori¬ 
entation  axis.  Its  positive  direction  is  upward.  Horizontal,  vertical,  up,  and  down  are 
relative  to  the  goniometer  base.  (It  is  assumed  that  the  Buerger  type  goniometer  is  used.) 

2.  The  orientation  axis  angle  9,  which  is  the  angle  between  the  positive  orientation 
axis  and  the  positive  horizontal  axis  lying  in  the  plane  of  the  .reciprocal  net.  The  positive 
horizontal  axis  is  directed  away  from  the  goniometer  head.  The  angle  is  expressed  in 
degrees. 

3.  The  precession  angle  ji  expressed  in  degrees. 


3.4  Category  III  -  Reciprocal  Net  biformation 
Card  8: 


columns 

1  -  3 

4-7 

8-10 

11  -  13 

14  -  25 

data 

LN 

NR 

LT 

PT 

SF 

format 

F3.0 

14 

F3.0 

F3.0 

E12.5 

Card  8  contains: 

1.  The  layer  number,  LN,  i.e.,  the  0th,  or  Ist,  or  2nd,  etc.,  net  above  the  zero 
level  net. 


2.  Number  (A  reflections  included  in  the  set,  NR. 
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3.  A  test  flag  LT  to  answer  the  question;  Is  this  the  last  set  of  data  in  this  setting 
to  be  read  into  the  computer? 

Yes:  LT  =  1 
No:  leave  blank. 

4.  A  test  flag  PT  to  indicate,  for  Weissenberg  data,  whether  the  compacted  or 
expanded  spots  were  read  In  this  set.  For  precession  or  integrated  Weissenberg  data, 
or  other  data  for  which  the  correction  does  not  apply,  the  space  is  left  blank. 

compacted:  PT  =  1 

ignore  correction:  leave  blank 
expanded:  PT  =  -  1. 

5.  If  two  standard  intensity  scales  have  been  used  to  estimate  intensities,  the  data 
estimated  with  one  of  these  scales  are  punched  as  negative  numbers.  This  field  contains 
the  factor  SF  by  which  the  negative  data  must  be  multiplied  to  bring  them  to  the  same 
scale  as  the  data  read  by  the  other  intensity  scale.  This  factor  is  punched  as  a  negative 
number.  If  only  one  standard  scale  has  been  used,  this  field  is  left  blank. 


Card  9: 


columns 

m 

4-6 

7-9 

10-12 

13-19 

20-26 

27-33 

34-40 

41-47 

48-54 

55-61 

data 

D 

a 

a 

€ 

6 

5 

4 

3 

2 

1 

AB 

format 

13 

13 

13 

13 

F7.0 

F7.0 

F7.0 

F7.0 

F7.0 

F7.0 

F7.5 

Card  9  contains  individual  reflection  data.  One  card  is  punched  for  each  reflection 
observed  in  the  set.  The  number  of  cards  in  this  set  must  equal  NR  as  specified  on  Card 
8  above.  The  card  contains  individual  reflection  data  as  follows; 

1.  The  Miller  indices  h,  k,  and  4. 

2.  The  integer  factor  c  which  is  the  factor  used  to  convert  E?  to  8?.  If  the  value 
of  €  is  1,  this  field  may  be  left  blank. 

3.  Six  fields  containing  the  intensity  readings  from  individual  films  either  from 
separate  exposures  or  from  one  or  more  film  packs.  There  can  be  a  maximum  of  six 
films.  The  fields  are  filled  from  ri^t  to  left  in  order  of  decreasing  exposure  of  the  films. 
The  sixth  field  from  the  left  must  contain  the  data  from  the  most  heavily  exposed  film. 

The  fifth  field  from  the  left  must  contain  the  data  from  the  secondmost  heavily  exposed 
film,  etc.  Where  intensities  are  too  strong  or  too  weak  to  be  estimated,  the  correspond¬ 
ing  field  is  left  blank,  even  when  preceded  or  followed  by  nonzero  fields.  Where  there 

are  fewer  than  six  films,  the  excess  fields  are  left  blank.  If  two  intensity  scales  are 
used  to  estimate  intensities,  the  data  estimated  from  one  of  the  scales  are  punched  as 
negative  numbers  (see  item  5  on  card  8  above). 

4.  An  absorption  correction  factor  AB.  This  is  the  factor  by  which  the  individual 
intensity  values  must  be  multiplied  to  correct  for  absorption.  If  there  are  to  be  no  cor¬ 
rections  for  absorption,  this  field  may  be  left  blank. 


14 
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3.5  Sequence  of  Cards 

Hie  data  deck  is  set  up  as  follows: 


General  Information: 

Card  1 

Card  2 

Card  3 

Card  4 

All  cards  5 

First  Setting: 

Card  6 

Card  7  (if  required) 

First  Set  of  Data: 

Card  8 

All  cards  9 

Second  Set  of  Data: 

etc. 

Card  8 

All  cards  9 

Second  Setting: 

Card  6 

Card  7  (if  required) 

First  Set  of  Data: 

etc. 

Card  8 

All  cards  9 

4.  OPERATING  PROCEDURE 

This  section  is  divided  into  two  parts.  Section  4.1  is  sufficient  (or  those  who  intend 
to  use  the  program  but  are  not  concerned  with  the  operating  procedure.  Section  4.2  will 
add  the  details  required  for  the  operation  of  the  program. 


4.1  Instructions  for  the  Crystallographer 

The  preceding  portions  of  this  exposition  should  be  carefully  studied  to  determine 
exactly  what  information  is  required,  and  how  it  is  to  he  presented  to  the  computer.  It 
would  be  wise  to  do  this  before  the  work  of  estimating  intensities  is  begun  so  that  the 
reflection  data  is  properly  indexed  and  is  recorded  in  a  form  convenient  for  key  punching. 
The  data  sheet  column  headings,  from  left  to  right,  should  be  h,  k,  c,  6,  5,  4,  3,  2,  1, 

AB  (cf.  formats  for  card  9  in  section  3).  The  absorption  correction  factors,  if  any,  should 
include  decimal  points.  If  no  absorption  correction  factors  are  to  be  applied,  their  fields 
may  be  left  blank.  If  Intensity  estimates  expressed  as  integer  numbers  are  ri^t  justified 
in  their  fields,  decimiil  points  need  not  be  punched.  Punched  decimal  points  may  appear 
anywhere  within  the  fields. 

Each  independently  observed  set  or  “layer”  of  data  must  be  handled  individually 
throughout  the  processing.  After  each  set  of  reflection  data  is  punched,  the  cards  are 
counted.  This  counting  should  be  done  mechanically  to  avoid  human  error. 

The  control  cards  (i.e.,  cards  other  than  reflection  cards)  may  be  inserted  in  the 
deck  after  the  individual  sets  of  reflection  cards  have  been  stacked  in  the  order  in  which 
they  are  to  be  loaded  into  the  computer.  The  information  on  the  control  cards  must  be 
consistent  with  this  ordering.  In  the  usual  case,  the  data  deck  will  be  loaded  onto  mag¬ 
netic  tape  for  input  to  the  computer. 
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It  is  suggested  that  the  data  processing  be  approached  in  the  following  way: 

Step  1.  The  data  processing  should  be  taken  through  (forced  through  if  necessary) 
completion  of  the  preliminary  output.  Any  gross  key  punching  or  transcribing  errors 
resulting  in  format  inconsistencies,  rejections  of  whole  sets  of  data,  or  abnormal  film 
factors  should  be  detected  at  this  stage.  The  preliminary  output  may  be  observed  on  the 
on-line  printer.  The  reflection  data  are  not  required  to  be  arranged  in  dictionary  order 
of  h,  k,  if  the  calculation  is  stopped  at  this  point.  Thus  the  first  part  of  the  calculation 
(DICNRY)  may  be  omitted. 

Step  2.  After  correction  for  the  errors  observed  above,  Step  1  is  repeated,  if  neces¬ 
sary,  in  order  to  find  all  remaining  detectable  errors. 

Step  3.  A  final  check  of  the  input  tape  (or  deck)  is  made,  by  processing  through  the 
output  of  rejected  reflections,  prior  to  committing  the  data  to  complete  processing.  Prior 
to  this  calculation,  it  would  be  wise  to  spread  the  film  factor  acceptance  limits  to  ensure 
that  no  additional  rejections  will  occur.  The  newly  computed  film  factors  will  be  differ¬ 
ent  from  the  old  ones  where  corrections  have  been  made.  Consequently,  some  data,  barely 
within  the  acceptable  range  in  the  prior  calculation,  may  lie  outside  of  that  range  in  the 
final  calculation. 

Step  4.  Finally,  the  processing  is  carried  through  to  the  farthest  stage  desired. 

The  BCD  output  tapes  should  be  set  aside  for  future  use.  If  corrections  to  these  tapes 
must  be  made,  it  will  te  necessary  to  have  the  data  on  them  punched  onto  cards.  The 
corrections  are  made  on  the  cards,  which  then  may  be  read  teck  onto  tape. 


4.2  Instructions  for  the  Machine  Operator 

The  data  reduction  is  divided  into  four  separate  calculations,  labeled  I,  II,  in,  and  IV, 
each  involving  Independent  read -in  of  its  own  program  deck.  These  four  calculations  are 
discussed  below  and  are  diagrammed  in  Appendix  A.  Each  successive  calculation  uses 
tapes  prepared  by  the  preceding  calculation.  To  gain  maximum  efficiency  the  four  cal¬ 
culations  should  be  done  in  immediate  succession. 

BCD  output  for  off-line  printing  is  on  tape  9.  For  XRDD,  this  output  also  may  be 
observed  on-line  by  depressing  sense  switch  6. 

Tables  1  and  2  present  summaries  of  the  tape  requirements  and  sense  switch  usage. 

The  complete  FORTRAN  source  program  listings  are  reproduced  in  Appendix  B. 


I.  DICNRY  (Dictionary  Ordering  of  Data) 

This  program  arranges  the  data  within  each  set  into  dictionary  order  of  h,  k,  and  -t. 
It  uses  tapes  4  and  5. 

Input  is  BCD  either  from  tape  4  or  from  cards  through  the  on-line  card  reader.  If 
the  on-line  card  reader  input  is  used,  depress  sense  switch  1. 

The  entire  program  output  is  on  tape  5.  This  output  is  mostly  binary  but  otherwise 
is  identical  with  the  input  tape  except  that  the  data  now  are  ordered  within  the  sets.  This 
tape  serves  as  input  to  the  next  calculation. 
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Table  1 

Summary  of  Tape  Requirements 


Program 

Tape"^ 

_ 

Input 

Output 

Inter- 

medigte 

Binary 

BCD 

I.  DICNRY 

4 

X 

X 

5 

X 

X 

X 

II.  XRDD 

4 

X 

X 

X 

5 

X 

X 

X 

5 

X 

X 

6 

X 

X 

X 

7 

X 

X 

7 

X 

X 

8 

X 

X 

X 

9t 

X 

X 

III.  SF  NORM 

4 

X 

X 

X 

5 

X 

X 

6 

X 

X 

7 

X 

X 

8 

X 

X 

97 

X 

X 

IV.  E  ORDER 

5 

X 

X 

9t 

X 

X 

’i'As  indicated  in  the  Table,  the  individual  tapes  may  serve  several 
widely  differing  functions  in  the  course  of  a  calculation.  Never¬ 
theless,  the  same  tapes  can  be  used  throughout.  In  practice,  how¬ 
ever,  the  user  may  wish  to  replace  and  save  tapes  4  (input),  5 
(input  and  output),  and  7  (output)  at  certain  stages  of  the  proc¬ 
essing  (for  example,  see  section  4.2,  Calculation  11). 
fTape  9  is  used  for  BCD  output  to  be  printed  off-line.  It  is  never 
rewound  by  the  program. 


Table  2 

Summary  of  Sense  Switch  Usage 


Sense  Switch 

Usage 

1 

Down  for  card  data  input  to  DICNRY  or  XRDD. 

2 

Down  for  BCD  tape  data  input  to  XRDD. 

3 

Down  if  data  processing  is  to  stop  after 
preliminary  output  is  completed. 

4 

Not  used. 

5 

Not  used. 

6 

Down  for  on-line  printing  of  output  from 

XRDD  (this  output  is  written  on  tape  9 
regardless  of  the  sense  switch  position). 
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The  dictionary  ordering  calculation  may  be  omitted  entirely  if  the  data  cards  are 
ordered  prior  to  loading  on  tape.  Ordinarily  It  will  be  omitted  In  the  preliminary  error 
detection  runs. 


II.  XRDD(Data  Reduction  through  the  Output  of  Structure  Factor  Data) 

This  program  uses  tapes  4,  5,  6,  7,  8,  and  9. 

Tape  5  is  used  for  input.  Tape  5  may  be  either  binary  or  BCD  depending  upon  whether 
or  not  Calculation  I  has  been  done.  If  it  is  BCD,  depress  sense  switch  2.  If  the  input  tape 
is  to  be  saved,  it  must  be  replaced,  or  dialed  out  and  another  tape  5  dialed  on,  after  the 
on-line  printer  reports  that  the  preliminary  averaging  and  Lp  corrections  have  been 
completed.  Usually,  no  interruption  in  the  calculation  is  required  for  this  replacement. 

Alternatively,  the  data  may  be  loaded  through  the  on-line  card  reader.  This  mode 
of  data  loading  is  selected  by  depressing  sense  switch  1. 

The  data  processing  proceeds  through  the  following  stages; 

1.  Read-in  of  each  set  of  data  into  the  computer,  calculation  of  average  intensities, 
application  of  Lp  corrections,  and  printout  of  film  factors.  This  printout  may  be  observed 
on-line  by  depressing  sense  switch  6. 

2.  Printout  of  rejected  data.  This  printout  may  be  observed  on-line  by  depressing 
sense  switch  6. 

Notification  of  the  completion  of  the  processing  to  this  stage  will  appear  on  the  on-line 
printer.  If  sense  switch  3  is  down,  calculation  will  stop  at  this  point.  To  continue  the 
calculation  press  the  START  button  on  the  console. 

3.  Least -squares  correlation  of  independently  observed  sets  of  data.  Notification  of 
completion  of  this  calculation  will  appear  on  the  on-line  printer. 

4.  Final  averaging  and  scaling.  Notification  of  completion  of  these  calculations  will 
appear  on  the  on-line  printer.  If  sense  switch  6  is  down  at  this  time,  the  on-line  printer 
will  print  notification  of  this  fact  and  the  calculation  will  stop.  To  continue,  with  sense 
switch  6  up  or  down,  push  START. 

5.  Structure  factor  output.  ITiis  may  be  observed  on-line  by  depressing  sense 
switch  6.  Notification  of  completion  of  the  output  will  appear  on  the  on-line  printer. 

Tapes  5  and  7  are  BCD  output  tapes  ready  for  input  to  the  Busing  least-squares 

refinement  program  (see  section  1.10).  Tape  5  contains  the  averaged  data  and  tape  7 
contains  the  unaveraged  data.  If  the  data  processing  is  to  be  extended  to  include  Calcu¬ 
lation  III  (calculation  and  output  of  normalized  structure  factors),  these  tapes  must  be 
replaced  and  saved  before  Calculation  111  is  begun. 

Tapes  4,  6,  and  8  are  used  as  input  to  Calculation  III. 


III.  SF  NORM  (Normalization  of  Structure  Factors) 

This  program  uses  tapes  4,  5,  6,  7,  8,  and  9.  Tapes  4,  6,  and  8  are  the  same  tapes  4, 
6,  and  8  involved  in  Calculation  II  above. 
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Tapes  5  and  7  are  BCD  output  tapes  ready  for  input  to  the  Fourier  program  of  Sly 
and  Shoemaker.  Tape  5  contains  the  averaged  data,  and  tape  7  contains  the  unaveraged 
data. 

Notification  that  both  the  calculation  and  output  have  been  completed  will  appear  on 
the  on-line  printer. 

Either  tape  S  or  7  may  be  used  as  input  to  Calculation  IV,  but  tape  5  usually  will  be 
used. 


IV.  E  ORDER  (Ordering  of  |EP-lData) 

This  program  uses  tape  5  for  input.  This  tape  is  either  tape  5  or  tape  7  from  the 
output  of  Calculation  111.  No  other  tapes  are  used.  The  contents  of  the  input  tape  are 
not  altered.  Output  is  on  tape  9. 

The  program  processes  5000  or  less  reflection  data  per  pass.  If  there  are  more 
than  5000  data  on  the  Input  tape,  the  program  processes  the  data  in  groups  of  5000  until 
all  are  done  (see  section  2.4). 

Notification  of  completion  of  the  calculation  will  appear  on  the  on-line  printer. 
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APPENDIX  A 


FLOW  DIAGRAMS  OF  THE  FOUR  SUCCESSIVE  CALCULATIONS 


The  four  diagrams  which  follow  indicate  the  main  stream  of  flow  of  the  data  process¬ 
ing  for  the  four  separate  calculations.  The  long-dash  boxes  enclose  the  portions  of  the 
calculations  under  control  of  the  FORTRAN  executive  programs.  FORTRAN  subroutines 
are  enclosed  in  short-dash  boxes.  The  subroutine  names  are  found  in  the  upper  left  cor¬ 
ners  of  the  boxes. 


Fig.  A1  •  Calculation  I 
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Sr  NORM 


READ  input 
FROM  TAPE  4 
(SCC  XROO  CKO 


calculate  and  output  normalized 
STRUCTURE  FACTORS  AND 
AVERAGES  FROM  AVERAGED  DATA 
IP  TAPE  BtFROM  ABSI) 

OP  TAPE  9.5 


place  unaveraged  data  on 
The  absolute  scale  calculate 
and  output  normalized  structure 
factors  for  the  UNAvERACED  DATA 
IP  TAPE  6(FR0M  LSOCOR) 

OP  TAPES 


Fig*  A3  -  Calculation  111 


Fig.  A4  -  Calculation  IV 


APPENDIX  B 


COMPLETE  FORTRAN  SOURCE  PROGRAM  LISTINGS 

The  complete  FORTRAN  source  program  listings  are  reproduced  in  this  appendix. 
FUNCTION  INTG  is  called  by  both  XRDD  EXC  and  SF  NORM. 
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C  niCNRY  FXC 

DIMENSION  0(6* 1000 1 f L 1 ( 1000 ) •L2 < 1000 1 *L3 (1000 ) iLPS ( 1000 ) tABSC 
l(1000)*AS(20) tASl(20) iRS(20) tBSl ( 20  I tCS ( 20 ) »AN (20 ) 

COMMON  LI •L2iL3tLPS*0*ABSC 
REWIND  A 
REWIND  3 

IFISENSE  SWITCH  11100»200 
100  READ  5 

READ  10iAR(BRtCRiALR»BER*GAR 
READ  20tRliR2 
READ  30tNASP«S21 

READ  40i(AS(I  )  lASK  n  tBSI  n*BSl(  I  )  >CS(  1  I  tANI  I  )  *I«1  tNASPI 
GO  TO  300 

200  READ  INPUT  TAPE  A«5 

READ  INPUT  TAPE  4i lOtARtBRtCRtALRtBERtGAR 
READ  INPUT  TAPE  4i20«Rl*R2 
READ  INPUT  TAPE  4« 30*NASP tS2 1 

READ  INPUT  TAPE  4<40» ( AS( I ) * ASl ( 1) tBS I  I ) • BSl ( 1 ) tCS( I ) • AN ( I  )  • I>1 * 
InASP) 

300  WRITE  OUTPUT  TAPE  9*5 

WRITE  OUTPUT  TAPE  5*10*AR*BR*CR*ALR»BER*GAR 
WRITE  TAPE  5  *R1*R2 

WRITE  TAPE  5  *NASP*S2I 

WRITE  TAPES  * ( AS( 1 ) tAS 1 ( I ) tBSI I ) tBSl ( I ) tCS ( I ) * AN ( I ) *  I *1  • 

iNASPt 

105  IF(SENSE  SWITCH  11110*210 
110  READ  S0*U1*U2*U3*CE*PE*REJ*YAM*S1NH 
IF(PF)910*310*120 
120  RFAD  A0*E1*F2*E3*PHI«UM 
GO  TO  310 

210  READ  INPUT  TAPE  4*50. U1 *U2.U3*CE *PE *REJ *  YAM *SI NH 
IF<PEt310*310*220 

220  READ  INPUT  TAPE  4*60*E1*E2*E3*PHI *UM 
310  WRITE  TAPE  5  *U1*"2*U3*CE*PE »REJ*YAM*SINH 

IF(PE1330*330*320 

320  WRITE  TAPE  5  *E 1 *E2 *E3*PH I *UM 

330  IF(SENSE  SWITCH  11130*230 
1'30  READ  70*SL.NR*TL.aSE*SFI 

READ  80*(L1 ( n*L2( I) *L3( I I.LPSI n  * IO( J*I ) * J* 1 *6 ) * ABSC ( I  1 *I>1*NR) 

GO  TO  340 

230  READ  INPUT  TAPE  4*70»SL*NR*TL*BSE*5FI 

READ  INPUT  TAPE  4*60* ( LI ( I ) *L2 ( 11 *L3 ( I ) *LPS( I ) * ( 01 J* I ) * J>1 *6 ) • 
lARSC(n*I«l*NR) 

340  WRITE  TAPE  5  *SL»NR»TL»8SE*SFI 

350  CALL  DICNRYINRl 
360  IF(TL)330*330*370 
370  IF(CE) 105*105*380 
380  REWIND  4 
REWIND  5 
PRINT  90 
PAUSE 

5  F0RMAT(72H 
1  1 
10  F0RMAT(6F11.7I 
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?0 

10 

40 

50 

60 

70 

80 

90 


FnRMAT(2E12.5l 
F0RMAT( I4»F5.4) 

F0RMAT( 5F8.4tF4.0) 

F0RMAT( 5F3.0.3F11.8  » 

F0RMAT(3F3.0«2F11.7> 

FORMAT ( F3« 0» I4»2F3.0*F12 •9) 

F0RMAT(413»6F7.0*F7.5» 

F0RMAT(62H  JOB  .FINISHEQ.,  IHP';T  FOR  DATA  REDUCTION  PROGRAM  IS  ON 
APF  5,  //////////  ) 

ENO(i«i»otO*n 
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SUBROUTINE  OICNRY(NR) 

DIMENSION  0(6*1000)  iLl  ( 1000 )  *1.2  ( 1000 )  iLS  ( 1000 )  iLPS  ( 1000 )  *ABSC 
1(1000) •M2( 1000) IMS) 1000) tKl ( 1000 ) tNl ( 1000 ) •N2 ( 1000 )• N3 ( 1000 ) *N) 
2  (1000)  *Q1  (6*1000)  tABSCK  1000)* J3(  100)  *1(2  (  100) 

COMMON  Ll*L2iL3*LPS*0*AaSC*Nl*N2*N3*NPSl*Ql*ABSCl 
FREQUENCY  100 ( 300 ) *  1 1 0 ( 2  *  1 *5 ) • 140( 300 ) t ISO ( 2  •  1 *6 ) *  180 ( 300 ) * 

1 190(5*1*3) *211 (0*1* S; *220(100 )*230(S*1 *2) *260(100) *270(3*1*2)* 
2300(100) *310 (5*1 iS) *331 (0*1*3) *340(23) *330(9*1 *2) *380(23) *390 
3 (5* 1*2) *430(5*1 *3) *470(23) *460(1*1*3) *470(1*1*9) *480(1*1*9) 

L»0 

MS1>L1 (  1  ) 

100  00  130  1-1*NR 

110  1F(L1(I)-MS1)120*130**30 

120  MSl-Ll ( I ) 

130  CONTINUE 
IGl-LKl) 

140  DO  170  Ial*NR 

130  1F(  IGl-LK  1  )  )160*170*  170 

160  IGl-LKl) 

170  CONTINUE 
175  MK-0 

IRO  no  210  I-ltNR 

190  1F( MS1-L1( I ) )210*200*210 

200  MK-MK4'1 

M2(MKI-L2( 1  ) 

M3  (MO -LSI  1  ) 

K1(MK)>1 

210  CONTINUE 

211  IF(  MO490*490*215 
215  MS2»M2(1) 

220  DO  250  I-1*MK 

230  1F(MS2-M2( I ) )250«23fl*240 

240  MS2-M2(I) 

250  CONTINUE 
IG2-M2 ( 1 ) 

260  DO  290  I>1»MK 

270  IF(M2(I)-lG2)290*290*2e0 

280  lG2aM2(I) 

290  CONTINUE 
295  ML-0 

300  DO  330  I-1*MK 

310  IF( MS2>M2( I ) )330*320*330 

320  KL-ML+1 

J3(ML)-M3( I ) 

K2(ML)-K1( 1 ) 

330  CONTINUE 

331  IF( ML)470»470*335 
335  MS3-J3(1) 

340  DO  370  I-1*ML 

350  1F(MS3«J3( I ) )370*370*360 

360  MS3-J3( 1 ) 

370  CONTINUE 
1G3-J3(1 ) 

380  DO  410  I-1*ML 
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390  1F( J3 ( I )“I63 )410»410»400 
400  IG3>J3( I ) 

410  CONTINUE 

420  DO  430  laltML 

430  IF( MS3-J3( 1 ) )4»0»440»450 

440  M>K2( I ) 

L»L+1 
N1(LI-MS1 
N2( L)>MS2 
N3( L)aHS3 
NPS1(L)-LPS(M) 

DO  44S  K>lt6 
443  01(K*L)aQ(IC«Mi 

ABSC1(L)»ABSC(M) 


430  CONTINUE 
HS3-HS341 

460  IF! 1G3-MS3)470»420»420 
470  MS2aMS24l 

480  1F< IG2-MS2 1490*293*293 
490  MS1-MS141 

300  IF( IGl-MSl »310*175*173 
310  WRITE  tape  3 

lABSCl(n*l*l*NR) 

80  F0RMAT(413  »6F7«1*F7« 
RETURN 

EN0(2«2«2*2*2I 


*(N1(  I  )  »N2(  n  *N3(  1  ).»NPS1  (l)*(Cil(<»I)»IL»l*6)» 
41 
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XROO  EXC 

DIMENSION  0(6*1000) *H1( 1000) •H2( 1000) *H3( 1000) t EPS ( lOOO) »AVEO( 1000 
1  )*REJO(1000).RENO(1000)*ABSC(1000).AVSE(S) *NPS(1000) *AS(20)*BS(20) 
2*CS(20) *AS1(20) •BSl(20)*AN(20)*SMyAM(10) •N1(1000)*N2(1000)* 

3  N3(1000) 

COMMON  H1*H2*H3«EPS*AVE0«REJ0*REN0*ABSC*0*N1*N2*N3*NPS 
FREQUENCY  lOS ( 0 *3* 1 ) » 123 (0*3 *1 ) • 153 (0*3* 1 ) *160(0*2 » 1 ) *172 ( 2 ) 

5  FORMAT (72H 
1  ) 

10  F0RMAT(6F11.7) 

20  FORMAT(2E12.3) 

30  FORMAT) I4*F5*4) 

40  FORMAT(3Fe.4*F4.0) 

30  F0RMAT(3F3.0i3F11.8) 

60  FORMAT(3F3*Ot2F11.7) 

70  FORMAT(F3*0*I4*2F3«0*E12*3) 

80  FORMAT(4F3.0*6F7.0*F7.5) 

500  FORMATdH  «41H  LEAST  SQUARES  LAYER  CORRELATION  FINISHED/////) 

510  FORMATdH  •25H  FINAL  AVERAGING  FINISHED/////) 

320  FORMATdH  *88H  DATA  REDUCTION  FINISHED*  RELEASE  SENSE  SWITCH  6  AND 
ICONTINUE  FOR  LOADING  OF  OUTPUT  TAPES/////) 

330  FORMATdHl) 

90  FORMATdH  •49HPREL I  MI  NARY  AVERAGING  AND  LP  CORRECTIONS  FINISHED/// 
1//) 

340  FORMAT dl3H  REPLACE  AND  SAVE  TAPES  3  AND  7.  READ  I 

IN  JOB  SF  NORM  FOR  CALC.  AND  OUTPUT  OF  E  AND  E**2-l *//////// 
2///) 

330  F0RHATdH0*3I4«8H  LAyERt3*I8) 

360  FORMAT (IHO) 

570  FORMATdHO*102H  5  4 

13  2  1) 

380  F0RMATdH0*18H  FILM  FACTORS  3E18.7) 

390  F0RMAT(68H  DATA  REDUCTION  FINISHED.  LOADING  OF  OUTPUT  TAPES  WILL 
INOW  PROCEED./////) 

REWIND  4 

IFISENSE  SWITCH  1)1000*900 

900  IFISENSE  SWITCH  2)902*901 
902  NNN«2 

GO  TO  1001 

1000  READ  5 

READ  10*ARtBR*CR*ALR*BER*GAR 
NNNa3 

GO  TO  1002 

901  NNN>1 

1001  REWIND  5 

READ  INPUT  TAPE  5*5 

READ  INPUT  TAPE  3*10*AR*BR*CR*ALR*BER*GAR 

1002  ALR»ALR»3. 14159/180.0 
BER«BFR»3. 14159/100.0 
GAR =GAR*3. 141 59/ 180.0 
COSAR=COSF(ALR I 
C0SBR»C0SF(RFR) 

COSGR«COSF(GAR) 

wR«AR*BR»CR*SORTFd.O-COSAR*COSAR-COSBR*LOSBR-COSGR*COSGR42**COSAR 
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1»C0SBR*C0SGR» 

SINAR»SINF(ALR) 

SINBR«SINF«BER) 

SINGR»SINF(GAR ) 

AA»BR*CR*SINAR/VR 

BB»AR*CR*SINRR/VR 

CC"AR*BR*SINGR/VR 

COSAL» ( COSRR*COSGR-COSAR ) / J  S 1 NBR*S I NGR ) 

COSBE*  (  COSAR*COSGR-COSRR  )  /  <  S  I  ^IAR*S  I  NGR  ) 

COSGA» ( COSAR*COSBR-COSGR ) / I S I NAR*S I NBR ) 

GO  TO  (903<1004»1003I«NNN 
1003  READ  20«R1«R2 

READ  30iNASPiS2l 

READ  40f ( AS( n  tASK I  ) tBSC 1) iBSl ( I ) .CS( I ) t AN( 1 ) *I>1 fNASP) 
GO  TO  1005 


1004  READ 

INPUT 

TAPE 

5*20tRl*R2 

READ 

INPUT 

TAPE 

5*30*NASP*S2I 

READ 

INPUT 

TAPE 

5*40t( AS( 1 ) tASl I n ,BS( I ) tBSl I  I ) «CS( I ) *AN( n  * 

1  I«l. 

NASP) 

GO  TO 

1005 

903  READ 

TAPE 

5  »R1»R2 

READ 

TAPE 

5  *NASP>S2I 

READ 

TAPE 

5  »( AS( I ) tASl ( n .BSI 1 ) tBSl ( I ) *CS( I ) tANI n  t 

1  I>ltNASP) 

1005  ZETA-0.0 
0NE23*0.0 
COSf,Al«0.0 
A«0«0 
B«0.0 

C  =  0,0 

A1«0,0 

A2»n,0 

A3=0.0 

NOZ»0 

EX»0«0 

SF»0.0 

N  1*0 

NSET»0 

NS  =  0 

BT1=0,0 

BT2«0.0 

NRTsD 

no  95  J=lt5 
95  AVSF(J)=0.0 
PRINT  530 
REWIND  7 

100  GO  TO(904»1007tl006)»NNN 

1006  READ  50»U1 tU2tU3*CEfPE*REJ*YAMt5INH 

GO  TO  1008 

1007  READ  INPUT  TAPE  5 t 5 0 . U 1 »U2 .U 3 » CE t PE • RE J i Y AM t S I NH 
GO  TO  1008 

904  READ  TAPE  5  tUl »U2 *U3 *CE »PE tRE J i Y AM »S I NH 

1008  NSET*NSFT-fl 
SMYA^'(NSET  )»YAM 
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105  IF( PC) 120»120»110 
no  60  T0(905il010fl009l  ,NNN 

1009  READ  60*E1  tE2*E3iPHI  ,ljM 
GO  TO  1011 

1010  READ  INPUT  TAPE  5 • 60 1 E 1 • E2 • E 3 •PH  1 »UM 
GO  TO  1011 

905  READ  TAPE  5  t E 1 . E 2 .E3 .PH I »UM 

1011  OHI=PHI*3.14159/lB0,0 
UM*UM*3,141 59/lHO.O 

120  GO  T0(906.1013»1012 ) .NNN 

1012  read  70..5L  tNRtTL.BSE.SFl 

READ  80. (HI ( I) tH2 ( 1 1 .Hai n .EPS(  1 ) . ( O ( J . 1  I . J  =  1 . 6  I . ABSC ( 1 )  .  I « 1 .  NR ) 
GO  TO  1014 

1013  READ  Input  tape  s.to.sl.nr.tl.bse.sfi 

READ  INPUT  TAPE  5 • 80 » ( Hi ( 1) .H2 (  I ) »H3 I  I ) . E PS (  I ) . ( 0 ( J • I ) » J  =  1 » 6 ) » 

1  ABSC( I liIaliNRI 
GO  TO  1014 

906  read  tape  5  .SL.NR.TL»BSE.SF1 

READ  TAPE  5  • ( Nl ( I  I • N2 (  11 tN 3 (I ) • NPS ( I ) • ( 0 IJ • 1) • J> 1 • 6 ) • 

1  ABSC ( 1) . I =1 .NR ) 

DO  907  1*1, NR 
Hl( I )=N1 ( I) 

H2( I )aN2( 1 1 
H3( I) ■N3{  1  ) 

IF(  APSCdl)  19.07,1906,1907 
1406  ABSC( II  *1.0 

1907  1F(NPS(II1  907.1908,907 

1908  NPS( I »  *1 

907  FPS( I )*NPS(  I  ) 

1014  MSLaiNTGI SL  1 
MUlslNTGdJl  I 
MU2»INTG(U2  ) 

MU3*1NTG(U3  ) 

NRT«NRT>MP 

WRITF  OUTPUT  TAPE  9,560 

WRITE  OUTPUT  TAPE  9 . 550 .HU  1 ,MU2 »MU3  ,MSL , NR 
IF(  <^PM<;F,5WtTCH  6)1015,1016 

1015  PRINT  560 

PRINT  550,MU1,MU2,MU3,MSL.NR 

1016  NS*NS+1 

call  XRIAVISFI ,NR,REJ,AVSE ) 

CALL  PSICAL (AA,Ufi,CC,COSAL.COS8E,COSGA,NR,YAM,Ul.U2,U3.SL»2ETA» 
10NE23,C0SGA1,A,8,C,A1 ,A2,A3) 

WRITE  OUTPUT  TAPE  9,570 

WRITE  OUTPUT  TAPE  9,  560 , ( AvSE ( J ) » J = 1  . 5 ) 

IF( SFMSFSWITCH  61  1017,125 

1017  PRINT  570 

PRINT  560.  (AVSFI J) .J=l  ,5) 

125  IF(PF)  140.140.no 

130  CALL  PRECLPIEI , E2 , E 3 , COSAL .COSBE ,C0SGA ,U 1 , U2 ,U3 ,PH I .COSGAl »NRf UM, 
1  YAM,0NE2  3,ZETA,A,B,C , A 1 , A2 , A3 , S I NH ) 

GO  TO  150 

140  call  WESSLP(YAM,NR.R1 ,R2,SINH,2ETA.6SE) 
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150 

155 

160 

170 


1701 

1702 


172 

175 

leo 

1  <50 


200 

210 

211 

212 


CALL  DIC0RD(U1 tU2*U3»CE»PE.SL»NR*TL ) 

IF( TL) 120tl20il60 
IF(CE)100.100il70 
CALL  REJECTINSl 
PRI^1T  90 

IF(SENSE  SWITCH  3I1701»1702 

PAUSF 

REWIND  5 

CALL  LSOCOR(NS) 

PRINT  500 

CALL  lAVOR (NStNOZ) 

PRINT  510 
SyAN  =  SMYAM(1  ) 

DO  180  I»1*NSET 

1F( SYAM-SMYAM( I ) ) 180  *  180 ♦ 1 75 

SYAMrSMYAM ( t ) 

CONTINUE 

call  ARSI  (SYAM.AS.flS.es*  ASltBSl  *NA.t,P,AN.S2  I  tNOZ.  ex*  SF.NI  .BT1.BT2I 

IF(SPNSF  SWITCH  61190.200 

PRINT  520 

PAUSF 

r,0  TO  210 

PRINT  590 

WRITF  OUTPUT  TAPE  9.530 
WRITF  OUTPUT  TAPE  9.5 
IF( SFNSFSWITCH  61211.212 
PRINT  550 
PRINT  5 

CALL  OUTPUT (NOZ.NS.N I .BT1.BT2.SF.EX.N0IB) 

WRITF  OUTPUT  TAPF  4.5 
WRITF  TAPF  4,NASP 

WRITE  TAPE  4.(AS( 1 ) .ASK  I) *dS( 1 ) .bSK 1) *CS( 1) *AN( 1 )*I«1*NASP) 

WRITE  TAPE  4.N02»NS.NI.BT1.BT2.SF.EX*N0IB*NRT 

REWIND  4 

PRINT  540 

PAUSF 

ENO(O.l.O.O.l) 
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SUBROUTINE  XRIAV  I SF 1 tNR tREJ (AVSE ) 

DIMENSION  Q(6«1000) iSE ( 9 t 1000 ) »AVSE I S ) (RE JO ( 1000 ) iRENO ( 1000 ) * 
IHKlOOOIt  AVEQ(IOOO)*  ABSC  ( 1000 )  tHa  (  1000  )  iHS  (1000  )  * 

1  EPS(IOOO) 

COMMON  Hl«H2fH3*EPStAVEO«REJOtRENO 
COMMON  ABSCtOiSE 

frequency  130(3001 *140 ( 1*4*12  I *60(1 »4 112) •200(300) *209 (0*1 tS I* 
1  20(0*1*31 *222 (300) *223 ( 1*0*20) *229(1*0*20) *249(300) * 

1  246(1*  6*  14)*289(300)*286(1*6*14) *319(1*6*14)* 

1  340(0*1*3)  *224(0*1*1)*329  (0*1*1) 

130  DO  199  1-1*NR 
DO  190  J  >1*9 

140  IF(Q(J*I))  170*190*180 
190  SE(J*I)  >  0*0 
160  60  TO  190 

170  Q(J*1)  >  SFI  •  0(J*I) 

ISO  SE(J*I)  ■  Q(J-^1*I)  /Q(J*I) 

80  IF(SE(J*1))  161  *  190*190 
181  SE( J*I )  >  SFI  •  SE( J*I ) 

190  CONTINUE 

90  IF(0(6*I ) )191*199*199 

191  0(6*1 )-0(6*I )*SF1 

199  CONTINUE 

00  270  J>1*9 

LUMP>0 

SIGMAl  >  0*0 
S1GMA2  >  0*0 

200  DO  220  1  ■  1*NR 

209  IF(SE(  J*n  )  220*220*210 

210  SIGMAl  >SE(J*I)  -^SIGMAl 
S1GMA2  >1*0  SIGMA2 

220  CONTINUE 

20  1F(SIGMA2)  21*21*221 

21  AVSE(J)«0*0 
60  TO  270 

221  AVSE(J)  >  SIGMAl  /  S1GMA2 

REJH  >  AVSE(J)  •  REJ 

REJL  >  AVSE(J)  /  REJ 

222  DO  240  I>1*NR 

223  IF(REJH  -SE(J*I))  230*230*224 

224  1F(SE( J*I )  )240*240*229 

229  IF(SE(J*I)  -  REJL  )  230*230*240 

230  SE(J*I)  >  -1*0 
LUMP>1 

240  CONTINUE 
IF(LUMP)270*270*241 

241  SIGMA1>0*0 

SIGMA2»0.0 

LUMP-0 

249  DO  260  I  >  1*NR 
246  1F(SE(J*I))  260*260*290 
290  SIGMAl  >  SE(J*1)  *  SIGMAl 

SIGI4A2  >  1*0  *  SICMA2 

260  CONTINUE 


32 


NAVAL  IISEARCH  LAIORATOR> 


GO  TO  221 
270  CONTINUE 

275  DO  280  J-l*4 
N  ■  J  ♦  1 

276  DO  280  L-  N*5 

280  AVSE(J)  ■  AVSE(J)  •  AVSE(t) 

285  DO  370  I  ■  ItNR 

SUHO  ■  0*0  t 

REJOin  ■  0*0 

RENOd)  >  0.0 

DO  310  J  *  1.5  4 

286  IF(SE(J«I))  287.290*300 

287  IF(J  -  1)  289*289.288 

288  1F(SE(J-1*1 >1289*289*300 

289  REJO(I)  ■  1.0 
GO  TO  310 

290  IF(Q(J*1 ) >310* 110*291 

291  IF( J-l>300*300*292 

292  1F(SE(J-1*I > >310*300*300 

300  SUMO  ■  0(J*I>  •  AVSE(J>  ♦  SUMO 

RENO(  I  >  ■  RENOd  >  ♦  1.0 
310  CONTINUE 

315  IF(SE(5*I>>  320*325*330 
320  REJOd>  ■  1*0 
GO  TO  340 

325  IF(0(6*I > >340*340*330 
330  SUMO  •  SUMO  •*  0(6*1  > 

RENOd  >  ■  RENOd  >  ♦  1.0 
340  lF(RENOd  >  >350*350*360 
350  AVEOd>  ■  0*0 
GO  TO  370 

360  AVEOd>  ■  SUMO  •ABSC(I>  /  RENOd  > 

370  CONTINUE 

FREOUENCY  267(0*l*6>*2e8(5*l*5>  *290( 0*10* 1 > *291  (0* 1  *6 > * 

1  292d*10*3)*90tl*4*12> 

RETURN 

ENO(1*1*0*0*0> 
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SUBROUTINE  PSICAL  ( AA *BB tCCtCOSAL *COSBE tCOSGAtNR* Y AM »U1 *U2 tUStSL* 
1  ZETA*ONE23»COS6Al*A*B*C«AltA2*A3i 

DIMENSION  XI ( 1000 ) •X2 ( 1000 )» X3 ( 1 000 1 (HI (1000 1 iH2( 1000) *83(1000 )• 
ISKIOOO)*  PSKIOOO)  *AVEO(1000)*REJO(1000)*RENO(1000)*Pt.(1000) 

1  «EPS(1000) 

COMMON  HI »H2  *H3  *EPS*AVEQ*REJQ»RENO 
COMMON  PSIiPLtSI 
COMMON  X1*X2*X3 

FREQUENCY  400 ( 1 *30* 10 ) *410 ( 1 *30* 10 ) •420( 1  * 30*10 1 *435 (300 ) « 

1  455(300)*47S(300)*495(300) 

A  •  AA  •  YAM 
B  »  BB  ♦  YAM 
C  •  CC  *  YAM 

All*  SQRTF( (U1»B)**2  ♦(U2*A)**2  ♦  2.0* (Ul*B )* 1 U2*A ) *COSGA ) 

A12*  SQRTF( (U3*B)*»2  ♦(U2*C)**2  ♦  2 .0* (U3*B »* ( U2*C ) •COSAL ) 

A21*  SQRTF( )U3*A)**2  ♦(U1*C)**2  ♦  2.0*IU3*A)»(U1*C)*C0SBE) 

A3  ■  SQRTF( (U1*A)**2  ♦(U2*B)**2  ♦(U3*C)**2 

1  2.0*(  (Ul*A)»(U2*B)*C0SGA  ♦  ( U1*A )  *  (U3*C)*C0S8E 

2  +(U2*BI*(U3*C)*C0SAL) ) 

BOOTCl  ■  Ul*U3*(A*A-C*C)*U2*U3*A»B*COSGA 
1  ♦(U3*U3  -  U1*U1)*A*C*  COSBE  -  U1*U2*B*C*  COSAL 

ADOTCl  »  U1*U2*(A*A  -  8*8  )  .  ■*•  (U2*U2  -  U1*U1 )  •A*B*COSGA 

1  -*■  U2*U3*A*C*COSBE  -  U1*U3*8*C*C0SAL 

A00TC2  «  U2*U3*(C*C  -  B*B  )  ♦  Ul*U2*A*t*C0SBE  -  U1*U3*A*B*C0S6A 

1  +(U2*U2  -  U3*U3)*8*C*COSAL 

UV  »  U1»U2 
U2W2  «  U1*U1  ♦  U3*U3 
VW  ■  U2*U3 
UW  ■  U1*U3 
U2V2  ■  Ul*Ul  ♦  U2*U2 
V2W2  ■  U2*U2  ♦  U3*U3 
UVW  «1.0/(Ul*Ul  ♦  U2*U2  ♦  U3*U3) 

400  1F(U1I430*410*430 
410  IF(U2)420*450*420 
420  1F(U3)430*470*430 
430  U0V,W«UVW/U1 

AOOTBl  «  VW  *A*A  -  UV  *A*C*COSBE  -  UW  *4*0*00304 
1  ♦U1*U1*B*C*C0SAL 
ONE23  •  -1.0 

COSALl  «  BOOTCl  /  (A21*  A3  ) 

COSBEl  «  AOOTCl  (  (All*  A3  ) 

COSGAl  «  AOOTBl  /  (All*  A21) 

A1  <  All 
A2  a  A21 

435  00  440  I  3  1»NR 


Xl(  n»(Hl(  n*UV  -  H2(l)*  U2W2  ♦  H3(I)*VW)*  UUVW 
X2( I )»(H1( I )*UW  +  M2(I|*VW  -  H3(I)*U2V2)*  UUVW 
440  X3( n»(Hl( I )*Ul  ♦  H2(I)*U2  ♦  H3(I)*U3)*  UVW 

GO  TO  490 

450  ADOTB2  =  -UV*C*C  +VW*A*C*COSBE  -  U3*U3* A*B*COSGA  +UW*B*C*C0SAL 
ONE23  ■  0.0 

COSALl  ■  BOOTCl  /(A21*A3  ) 

COSBEl  «  AOOTC2  /(A12*A3  ) 

COSGAl  ■  AOOTB2  /(A12*A21) 
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WUVW>UVW/U3 
A1  ■  A12 
A2  ■  A21 

ASS  D0A60  I>  1*NR 

X1(I)«  IHl(n»UV  -  H2(I)*U2W2  ♦  H3(I)*VW|*  WUVW 
X2(n>  (H1(I)«V2W2  -  H2(I)«UV  -  H3(n*UW)*  WUVW 
A60  X3(n>  U1  *  H2(II»U2  *  H3(I)*U3)*  UVW 

GO  TO  490 

470  AOOTB3  •  -UW*B*8  ♦  VW«A*B*C0S0A  -  U2*U2*A*C*COSBE  ♦UV*B«C*COSAL 

ONE23  >  ItO 

COSALl  >  -  A00TC2/(A12*  A3) 

COSBEl  >  A00TCI/(A11«  A3) 

COSGAl  >  AOOTB3/(A11«  A12 ) 

VUVW>UVW/U2 
A1  >  All 
A2  >  A12 

473  00480  I  •  1*NR 

XKI)  *  (H1(I)»V2W2  -  H2(1)*UV  -  H3I1)*UW)*  VUVW 

X2(I)  •  (H1(1)*UW  H2(I)»VW  -  H3I1)*U2V2)*  VUVW 
480  X3(l)  >  (H1<I)»U1  *  H2(I)»U2  *  HStI)»US)*  UVW 
490  SIN2A1  >  1.0  *  C0SAL1»C0SAL1 

SIN2B1  >  1.0  •  C0SBE1«C0SBE1 

SIN2G1  >  1.0  -  COSGAl •COSGAl 

ZETA  >(SL*A3)**2  •  SIN2B1*(1*0  -(  COSALl  -  C0SBE1*C0SGA1 )**2 

1  /(SIN2G1*S1N2B1>) 

DEL  •  SL*A3*(StN2Gl«COSaEl  *  C0S6A1*(C0SAL1  •  C0SGA1*C0SBE1 ) t / 

1  S1N2G1 

EPSI  ■  SL*A3*(C0SAL1-C0SGA1«C0SBE1)/SIN2G1 

499  00  900  1-1 *NR 

IP(X3( I ) )496«900*900 

496  XKD— XKI) 

X2(I)— X2(I) 

X3(  I  )><-X3(  1 ) 

500  PSKl)  •(XKI)*A1  ♦  OEL)»*2  •  (X2tl)»A2  ♦  EPSI)»*2 

1  •  2«0*(X1(1)*A1-»OEL)«(X2II)«A2  •«>  EPSI)*C0S6A1 

RETURN 

ENOI1*1*0*0»0) 
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SUBROUTINE  PRECUPIEl*E2tES*  C0SAL»C0SBE»C0S6A*U1*U2»U3* 

1  PHl»COS6AltNR*UHtYAH*ONE23*2ETA*A*BtC«Al*A2«A3»SINH)  | 

DIMENSION  COST(1000)(SlNT(1000)*XltlOOO)*X2(1000)*X3(1000>t 

1  PSII1000)»PL(1000)«SI(1000)*AVEO(1000)*Hl(1000)tH2(1000)f  l 

1  H3(1000)«EPS(1000)*REJO(1000)*RENOI1000)  I 

FREQUENCY  5 ( 2*1 *1 ) • 30( 1 10» 1 ) *6S( 2*1  *  1) *80( 1 *0*1 ) »90< 1 tOf 1 ) t 

1  145(30C)*163(2*1*1)*180<1*0*1I»190(1»0»1)»2S»(300)(  I 

1  273(300) 

COMMON  Hl*H2*H3*EPS«AVEQtREJO*RENO  | 

COMMON  PSItPLtSl 
COMMON  Xl»X2*X3 
COMMON  COST fS I NT 

VER  ■SQRTF((E1»A)«(EI»A)  ♦  (E2*B)*(E2*B)  ♦  (E3*C)*(E3*C) 

1  •*  2.0*( (El*AI«(E2«B)»COS6A  ♦(E1»A)*(E3*C  )*C0SBE  * 

1  (E2»B)*(E3*C)«COSAL)) 

VOOTAl»El«U2*A*A  -  E2«U1*B»B  ♦  (E2*U2  -  E1«U1 )*A«B*COSGA 
1  ♦  E3*U2»A*C»COSBE  -  C3«U1«B*C*C0SAL 

VOOTBl*  E1»U3*A*A  -  E3*U1»C*C  ♦  E2»U3*A*BftCOS6A 
1  ♦(E3*U3  -  E1*U1)»A*C*C0SBE  •  E2*U1«8»C»C0SAL 

VDOTA2>  -  E2*U3*B*8  *  E3*U2*C*C  -  E1*U3*A«B»C0S6A 
1  •*  (£2*02  •  E3«U3)*B*C«COSAL  *  E 1*U2»A*C*C0S8E 

3  1F(0NE23)  10*20*30 
10  VOOTA  >  VOOTAl 
VOOTB  ■  VOOTBl 
GO  TO  40 

20  VOOTA  >  VD0TA2 
VOOTB  >  VOOTBl 
GO  TO  40 

30  VOOTA  >  VOOTAl 
VOOTB  ■  -  VOOTA2 
40  1F(1.0«ABSF(COSGA1) )42*42«43 

42  StNGAl«0*0 
GO  TO  46 

43  SlNGAl  ■  SQRTF  <1*0  •  COSGAl*  COSGAl) 

46  COSVA  >V00TA/(A1*  VER) 

COSVB  >V00TB/(A2  *  VER) 

CVAl  »  E3*(U1*U1  ♦U2*U2)  -  U3«(U1*E1  ♦  U2»E2) 

CVBl  ■  U2*(U3*E3  *  U1*E1)  •  E2*(U1»U1  ♦  U3*U3) 

CVA2  •  U1*(UI»E2  ♦  U3*E3)  -  E1*(U2*U2  ♦  U3»U3) 

COSPH  ■  COSF(PHl) 

SINPH  >  SINF(PHI) 

30  IF (ABSF( COSVA)  •  ABSF(COSVB))  60*160*160 
60  IF (1.0-ABSF( COSVB) >62*62*63 

62  SINVB>0*0 
GO  TO  66 

63  SINVB  >  SQRTF (  1*0  •  COSVB  •  COSVB) 

66  CAPHVB  ■  COSPH  *  COSVB  •  SINPH  •  SINVB 
SAPHVB  ■  SINPH  *  COSVB  ♦  COSPH  ♦  SINVB 
CDPHVB  »  COSPH  *  COSVB  +  SINPH  *  SINVB 
SOPHVB  ■  SINPH  •  COSVB  -  COSPH  *  SINVB 
63  IF(ONE23)  70*  70*92 
70  IF(CVBl)  90  *  80*80 
80  IF(COSVB)110*100*100 
90  IF(COSVB) 130*120*120 
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92  1F(CVA2)  80*80*90 

100  COSA  ■  CAPHVB  •  C0S6A1  ♦  SAPHV6*  S1N6A1 

SINA  •  >  SAPHVB  •  C0S6A1  •*  CAPHVB*  SINGAl 

60  TO  140 

110  COSA  »CAPHVB  »  C0S6A1  >  SAPHVB  •  SINGAl 
SINB  ■  SAPHVB  «  C0S6A1  •  CAPHVB  *  SINGAl 
60  TO  140 

120  COSA  -  CDPHVB  *  COSGAl  •*  SDPHVB  *  SINGAl 

SINA  >  -SOPHVB  •  COSGAl  ♦  CDPHVB  *  SINGAl 

GO  TO  140 

130  COSA  >-C0PHVB  •  COSGAl  -  SDPHVB  •  SINGAl 
SINA  ■  SDPHVB  «  COSGAl  >  CDPHVB  •  SINGAl 
140  CAT  >  A2  •  SINGAl  *  SINA 
CBT  ■  C0SA/(2*0  •  All 

SAT  ■  A2  •  SINGAl*  COSA 
SBT  ■  SINA/(2.0  *  Al) 

A12  >  A1*A1 

A22  >  A2*A2 

149  00  ISO  I  >  1*NP 

XY2  >  IX2(n*X2(II*A22  -  XI  ( li  *X1 1 1 )  *A12  •  PSKIII/  XKII 

COST! I)  ■  X2(II  *CAT  -  XY2  *  CBT 
190  SINTin  >  X2(n  *  SAT  4  XYZ  *  SBT 
60  TO  270 

160  IFd.O-ABSFtCOSVAtl  162*162*163 

162  SINVA>0*0 
60  TO  164 

163  SINVA  ■  SORTF(  1.0  •  COSVA*COSVA) 

164  CAPHVA  >  COSPH  *  COSVA  -  SINPH*  SINVA 

SAPHVA  «  SINPH  *  COSVA  •*  COSPH  *  SINVA 

COPHVA  >  COSPH  *  COSVA  4  SINPH  *  SINVA 

SOPHVA  >  SINPH  *  COSVA  •  COSPH  *  SINVA 

169  1F(0NE23>  170*200*170 

170  IFICVAll  190*180*180 
180  IFICOSVAI  220*210*210 
190  IFICOSVAi  240*230*230 
200  IF(CVA2I  190*180*180 

210  COSA  >  CAPHVA  *  COSGAl  -  SAPHVA  *  SINGAl 

SINA  ■  SAPHVA  *  COSGAl  4  CAPHVA  *  SINGAl 

60  TO  290 

220  COSA  -  -  CAPHVA  *  COSGAl  4  SAPHVA  *  SINGAl 

SINA  ■  -  SAPHVA  •  COSGAl  -  CAPHVA  •  SINGAl 

60  TO  290 

230  COSA  «  COPHVA  *  COSGAl  -  SDPHVA  *  SINGAl 

SINA  •  SOPHVA  *  COSGAl  4  COPHVA  *  SINGAl 

GO  TO  250 

240  COSA  ■  -  COPHVA  *  COSGAl  4  SOPHVA  *  SINGAl 

SINA  «  -  SDPHVA  *  COSGAl  -  COPHVA  *  SINGAl 

290  CAT  «  A1*SINGA1*SINA 
CBT  ■  C0SA/(2.0*  A2> 

SAT  »  Al  *  SINGAl  •  COSA 
SBT  ■  SlNA/<  2.0*  A2I 
A12  -  A1*A1 
A22  ■  A2*A2 
299  DO  260  I  >  1*  NR 
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XYZ  -  (xin)*xi(n  *  A12  -  X2(n*  X2(i)  *  A22  -  pskd)/  x2(I) 
COST (11  ■  XKI)  *  CAT  -  XYZ  *  CBT 
260  SINTCII  ■  XKI)  •  SAT  ♦  XYZ  *  SBT 

270  SINU  ■  SINF(UM) 

COSU  ■  COSF(UM) 

TAN2U  ■(SINU/C0SU)»»2 
SUCU  >  SINU*COSU 
IF(ZETA)271t272*272 

271  ZETA-OiO 

272  ETAA  ■  SINU*SINU  -  1*0  ♦(COSU  -SQRTF(ZETA( )**2 

ETAB  ■  2«0»  SINU 

YAMB  ■  (1.0/(2.0  •  YAM))*«2 


275  DO  280  I  >IiNR 

COSETA  ■  (PSI ( n  ♦  ETAA)/ 

C02ETA  ■  COSETA  *  COSETA 

SI2ETA  -PSKI)  -C02ETA 

IF(SI2ETA)276»276»277 

276  SINETA>0.0 

REJQ( I )>REJO( I)  ♦  10«0 


ETAB 


GO  TO  278 

277  SINETA  ■  SORTF(SI2ETA) 

278  ABC  ■  C02ETA  *  SINT(II«  SINT(I)  ♦  SI2ETA  •COST(I)»  COST(I) 

OEF  ■  2.0*  SINETA*COSETA  *  SINTd)*  COST(I) 

PL(I)«  SINETA  •  SUCU  /(4.0  ♦(  2.0  -  (PSI(I)+  2ETA)*(PSI ( I )*ZETA) ) 
1  )/  (  1.0/  (1.0  ♦  TAN2U*IABC  ♦  DEF  ) /( PSI ( 1 )*PS1 ( I ) ) ) 

1  *1.0/  (1.0  ♦  TAN2U*(AeC  -  DEF  ) / (PSI 1 1 )*PS1 1 1)1 1 ) 

1  *50. 

AVEO( 1 )  >  AVEQ( 1 )  *  PL(  1 ) 

PL(I)  -  SQRTF(PL(II) 

Sl(n«(PSI  (n+ZETA)*YAMB 
1F( SlNH'Sl ( I ) 1279*279 *280 

279  AVEQ( 1 )*AVEQ( I )*1.9 

280  CONTINUE 
RETURN 

eN0(l*l*0*0*0) 
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SUBROUTINE  WESSLP<  YAM»NR*Rl«R2*SlNHt  ZETA»BSC) 

DIMENSION  AVEQ< 1000 ) *PS1 I 1000 ) *PL(1000I *it ( 1000) (HI I 1000) • 

1  H2 ( 1000 ) tHS ( 1000 ) (EPSl 1000 ) *REJQ( 1000) (RENOI 1000) 

COMMON  Hl*H2tH3»EPS*AVEO»REJO*RENO 
COMMON  PSItPLtSI 

FREQUENCY  90(250 )• 10510*4*1 ) *120( 2 *2 tl ) •1901290) •150(250) 

1 *91 (1*1 *100) *139 (100*10*1) *155 (100* 10^1 )^U1(100^1^1) 

YAMB  >(l*0/(2*0  •  YAM))««2 

90  DO  100  I  -  1*NR 

91  IF(4.0-  P$1(!)-ZETA)  92*92*93 

92  REJO(I)  ■  REJO(l)  ♦  10*0 
PL(1)  >  0*0 

60  TO  94 

93  PL(1)  ■  2«0*  SORTF(PSI(I)«I4*0-PS1(1)  -  ZETA))/ 

1  (4*0  ••  (2«0  -  PSI(t)  -  ZETA)»(2*0  -  PSl(l)  -  ZETA)) 

94  AVEO(I)  >  AVEQ(l)  •  PL(1) 

PL( I )  >  SORTF(PL( 1) ) 

SKI)  >  (PSKI)  *  ZETA)  •  YAMB 

1F(SINH-SI(1))95*95*100 
99  AVEO(l)aAVEO(I)*l*9 
100  CONTINUE 

109  1F(ZETA>1«OE-10)1TO*170*110 

110  PSIM  -  4*0  >  ZETA 

RZ  >  100*0  •  SQRTF( ZETA)/ (12 *5*036  *(2*0  •  R1/S0RTF(PS1M) 

120  1F(SSE)150*170*130 
130  DO  1451  •  1*NR 

139  1F(REJQ(I)-10*0)140*149*149 

140  YDEL  •  RZ  •  SORTF(  PSIM/PSK  1 )  -  1*0) 

141  IF(l*0-VOEI.-l*0E-02)142*142*143 

142  AVEQ(l)a0*0 
REJO(1)>REJO(I)«20*0 
60  TO  149 

143  AVEO(l)  >  AVEQ(l)  •  (1«0-  YOCt) 

PL  ( 1 )  "PL  ( 1 )  •SORTF  ( 1  *0-Y0f  I. ) 

145  CONTINUE 
60  TO  170 

150  DO  169  I  •  1*NR 

195  IF(REJQ(1)-10*0)160*165*165 

160  YDEL  >  RZ  •  S0RTF(PSIN/PS1 ( 1 )  -1*0) 

AVEQ(I)  -  AVEQ(I)  •  (1*0  ♦  YDEL) 

PL( I )>PL( I )*SQRTF(1*O^YOEL) 

169  CONTINUE 

170  RETURN 
ENO(2*2*2*2*2) 
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SUBROUTINE  DICORD(Ul»U2*U3*CEtPE»SL*NR*TL) 

DIMENSION  HI ( 1000) »H2 ( 1000 ) fH3( 1000) tEPSI 1000) tAVEOI 1000) * 
1REJQ( 1000 ) tRENOI 1000 ) *PL ( 1000 ) tSI 11000) tPSl ( 1000) 

COMMON  HI iH2*H3  tEPS«AVEQ*REJO(RENO 
COMMON  PSIiPL»SI 
WRITE  TAPE  7»Ul*U2*U3tCEfPE 
WRITE  TAPE  7»SLiNR»TL 

WRITE  TAPE  7*(Hl(I)tH2(I)«H3(I)*EPSn)tAVEO(I)«REJO(I)tRENO(I)* 
lPL(l)*SI(l)«I>ltNR) 

RETURN 

EN0(0*1*0*0*0) 
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SUBROUTINE  REJECT (NS) 

DIMENSION  Hi ( 10001 •H2 ( 1000 ) *H3 ( 1000 ) •AVE0( 1000) *KE JO (1000) *PL( 1000 
1 ) fSI ( 1000) »SPE( 1000) *PSE( 1000 )• ESP ( 1000) 

COMMON  H1*H2(H3*SPE»AVEQ*KEJ0*PSE»ESP*PL>SI 
FREQUENCY  3 ( 20 ) »10( 230 ) *20(0*20 » 1 ) 

REWIND  7 

IFCSENSESWITCM  6)2*3 

2  PRINT  500 
PRINT  560 

3  WRITE  OUTPUT  TAPE  9*500 
WRITE  OUTPUT  TAPE  9*  560 

5  DO  50  J-1*NS 

READ  TAPE  7*U1 *U2*U3*CE*P£ 

READ  TAPE  7*SL*NR*TL 

READ  TAPE  7  * ( HI ( I ) •H2 ( I ) •H3( 1 ) *EPS*AVE0( I ) *RE JO( 1 ) tEPStPC ( I ) * 

1  SI(I)*I-1*NR) 

MSL>1NTC(SU) 

MU1>INT6(U1 ) 

MU2>INT6(U2 ) 

MU3«INT6(U3) 

IF(SENSESW1TCH  6)6*7 

6  PRINT  560 

PRINT  670*MU1*MU2*MU3*MSL 

7  WRITE  OUTPUT  TAPE  9*560 

WRITE  OUTPUT  TAPE  9*  670»MU1*)>(U2»MU3»MSL 
10  00  40  I>1*NR 
20  IFIREJQd))  40*40*30 
30  MH1*INTG(H1 (1  )) 

MH2*1NTG(H2(I)) 

MH3*|NTG(H3(I)) 

MRJaINTG(REJO( I ) ) 

IF(SENSESWITCH  6)32*33 

32  PRINT  510*MHl*MH2*MH3*SI(I)*Pt.lI)*AVEO(I)*MRj 

33  WRITE  OUTPUT  TAPE  9*  510*HH1 »MH2 »MH3»SI ( I ) *Pl ( I ) *AVEO( I ) *MRJ 
40  CONTINUE 

50  CONTINUE 
REWIND  7 
RETURN 

500  FORMAT (1H1*84H  H  K  I.  SIN()/L  LP  FACTOR 

1  F4*2  R) 

560  FORMAT! IHO) 

670  FORMAT! 1H0*3I3*9H  LAYCRI3) 

510  F0RMAT!1H0*2H  3I4«F19*7*F18*6»F19.3»I 14) 

EN0(2«2*2*2*2) 
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SUBROUTINE  LSOCORINS) 

DIMENSION  HI ( 1000) •H2( 1000) tHSt 1000) tREJOl 1000) *AVEO( 1000) t 
161 ( 1000 ) *62 ( 1000 ) *63 ( 1000 ) tAVEOl ( 1000 ) «RE JOl 1 1000 ) vRENOl 1 1000 ) • 
2CI30)tA(30*50)»M4ll00)*M5llO0)*MS(100)»TRENOI100)*TAVEOa0O)* 

SN4(  100 )»NS 1 100  )*NS(  100)  tTAVEOK  100)  tTRENOK  100)  tXtSO  )• 

4REN0( 1000 ) tCPS) 1000 ) tPL ( 1000 ) *SI 1 1000 ) *B I  90 ) 

COMMON  HltH2*H3»EPS*AVEQ*REJO*RENO 
COMMON  M4»M9*X«C(TREN0*TAVE0*N4'*N5*N6»M6*TAVE01 
COMMON  PL*S1 *61(62 •63*AVEQ1«REJ01»REN01»A*B»TREN01 
FREQUENCY  270(20) »292(0«1*1)*299(0*6*1) »320(30*1*30) * 

1330(0*1 *10) *331 (0*9*1) *319 (300) *399 (30*1 *30) *360(0*1 *10)* 

2361 (0*9*1) *390(300) *379 (1*2*1) *410(1 *19*1) *4201 1*9*1) *4301 1*9*1)* 
3460(1 *1*4 >*470(1*1 *4) *471(1*1*4) *472(1 *1*4) *491 (0*6*1) •499(6*2*11 • 
4926(300) 

REWIND  6 
REWIND  7 
REWIND  8 
DO  200  1-1 *NS 

200  A(1*1)>0«0 
JK-1 

Nl-3*NS-3 

READ  TAPE  7*U1*U2*U3*CE*PE 
IF(CE)  20*20*201 

201  DO  202  1>1*NS 

202  X(l)  -IaO 
REWIND  7 
60  TO  929 

269  JK  >  JK  *  1 

270  DO  260  !•  1*N1 

2S0  BACK  SPACE  7 

N1  ■  N1  -  3 

290  READ  TAPE  7*  U1*U2*U3*CE*PE 
20  READ  TAPE  7*  SL*NR*TL 

READ  TAPE  7*(H1 ( 1 ) tH2( 1 ) *N3( I ) *  PS*  AVEO( I ) *REJO( 1 ) *RENO( I ) * 

1  PS*  PS*I  •  1*NR) 

292  1F(CC)291*291*900 

291  N  * 

299  1F(TL)  300*300*310 

300  READ  TAPE  7*  V1*V2*V3*CE1*PE1 
READ  TAPE  7*  SL1*NR1*T(. 

READ  TAPE  7* (61( 1 ) *62 ( 1 ) *63( 1 ) *  PS*  AVEOl ( 1 ) *REJ01 ( I ) (RENOK 1 ) * 

1  PS*  PS*I  •1*NR1) 

N  ■  N  *1 
A(JK*N)>0*0 
A(N*JK)>0*0 
60  TO  299 

310  READ  TAPE  7*  V1*V2*V3*CE1*PE1 
READ  TAPE  7*  SL1*NR1*TL 

READ  TAPE  7  * ( 6l ( 1 ) *62 ( 1 1 *63( 1 ) *  PS*  AVEOl ( I ) *REJ01 ( I ) *REN01 ( I ) * 

1  PS*  PS*I  -l*Nf(l) 

N  ■  N  ♦  1 

1F(U1-V1)314*311*314 

311  1F(U2-V2)914*912*314 

312  1F(U3>V3)314*313*314 
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313  A(JKtN)>0«0 
A(N»JK)>0«0 
GO  TO  491 

314  SRENO  -  0*0 
SAVEQ  >  0*0 
SAVEQl*  0*0 
J1  >  0 

319  DO  340  I  ■  1*NR 

320  IF(ABSP.(Si.l-ABSF(Vl»Hlcn4V2*H2(Il4V3*H3(llM-l«0C<-09}3i0«3»0*34e 

330  IFiAVEQd)  )331«331«332 

331  IFfREJOd)  I332*332«340 

332  J1  -  J1  ■«>  1 

M4(  Jl)>INTG(Hldn 
M9(JliaINTG(H2( 1)1 
M6(Jl)>lNTG(H3dtl 
TRENO(Jl)  >  RENOdI 
TAVEO(Jl)  •  AVEQdt 
340  CONTINUE 
J2  >  0 

390  00  370  1  •  ItNRl 

399  IFUBSFfSt  -ABSFCUl*Gld  l♦U2•G2 d  KU3*G3d  1 1  )-1.0EH)9)3G0*3A0»370 

360  IFlAVEOld  )1361*361»362 

361  IFlREJQld  ))362«362*3T0 

362  J2  >  J2  ♦  1 
N4IJ21alNTG(61(lll 
N9(J2)>lNTG(G2dll 
N6U2l«lNTG(63dll 
TREN01(J2)>  RENOld) 

TAVE01(J2|a  AVEOldl 

370  CONTINUE 
L  ■  0 
M  ■  0 

379  IFIJI  -  J21  390t3«0*3B0 
390  J  >  J2 
JBaJl 
GO  TO  400 
390  J  >  J1 
JB>J2 

400  00  469  I  ■  1*J 

409  1L>I->L 
1N"14N 

410  IFIM4dL)-N4dNll440*420»490 
420  XF(M9dl.l<-N9dNI)440»43e*690 
430  lF(M6dL)-N6dMn440*471»490 
440  L  -  L  ♦  1 

GO  TO  460 
490  N  >  N  ♦  1 

460  1F(JB-1L  )490«470*470 

470  IF(JB-1M  )490*409*409 

471  IF< JB-IL)490»472»472 

472  IF( JB-IM)490«4a0f460 

460  SRENO  •  SRENO  *  TRENOdt  1  ♦  TRENOldN  I 

SAVEQ  •  SAVEO  *  (TRENOdL  I  4  TRENOldN  !)•  TAVEQdL  1/ 

KTAVEOdL  I  *  TAVEOldM  1) 
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ABB  SAVEQl  >  SAVEOl  ♦  (TRENOIll.  )  •*’  TRENOKIM  )|*  TAVEOlIlM  »/ 
KTAVEQCIL  >  ♦  TAVEQIIIM  II 

490  A(JK*JKI  >  AIJKtJKI  ♦  SAVEQ  •  SAVEO  /  BRENO 

AlNtNI  >  AINtNI  SAVCQl  •  BAVEOl  /  SRENO 

Af JKtNI  >  -  SAVEQ  •  SAVCQl  /  SRCNO 
AINvJKI  ■  A(JR*NI 

491  IFITL)  BlOtSlOt  499 
499  IPICEllBlOfSlO*  8S9 

900  REWIND  7 
SN-NS 
AVA-0.0 

DO  901  1«1»NS 

901  AVAaAVA4A(l«l I 
AVA>AVA/SN 

DO  902  ial«NS 
DO  902  J>1*NS 

902  ACl*J|aA(I»JI/AVA 

WRITE  TAPE  S*<(An«JI«l>l*NSI«J-l*NSI 
REWIND  4 

CALL  EI6EN(A*C*NSI 
SMCaCIll 
DO  919  la2*NS 
IP(C(1I-SMC 1912 *919 *919 
912  SMCaCdl 
919  CONTINUE 

READ  TAPE  4*( lAI 1«JI • laltNSI •J*1»NSI 
REWIND  4 

CALL  VCTR(A«X*NS*SNC) 

929  DO  940  1  ■  1*NS 

READ  TAPE  7*  V1*V2»V9 •CEltPEl 
READ  TAPE  7*  SL  (NR  tTL 

READ  TAPE  T*IH1(JI«H2(4I«H9( JI*CPSIJ|»AVCOU>*REJO( JI»RENOUl* 
lPL(J|«Sl<U|*J>lfNR| 

924  DO  990  J  •  1*NR 
990  AVCOtJI  ■  X(II  •  AVCOtJI 
WRITE  TAPE  StVl»V2»V9«CEl*PCl 
WRITE  TAPE  4*SL»NR»TL 

WRITE  tape  4t(HlUI*M2(JI*H9CJI*EPSIJ|»AVEO(J|»REJOIJItRENOUI* 
lPLIJ|»SI(JI»Jal»NR| 

WRITE  TAPE  BfVl»V2fV9«CEl*PEl 
WRITE  TAPE  B*SLfNR»TL 

940  WRITE  TAPE  B*IH1( J| •H2(UI (HBIJI tEPSI Jl •AVEOl Jl •RCJOt Jl •RENOI Jl • 
1PL(  JIvSIUlvJalcNRl 
949  REWIND  4 
REWIND  7 
REWIND  B 
RETURN 

CNOIO*1*««0»OI 
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SUBROUTINE  I  AVON  ( NS  »NOZ )  ^ 

DIMENSION  HI (1000 t»H2( 1000) *>13(1 000) »ERS( 1000) »AVE0( 1000) (REJOI 100 
1  0)  vRENO ( 1000)  tSK  1000)  tHU)  50)  tHlHCSO)  tEKlOOO)  (FE ( 1000 )  tFS  ( 1000) 
1  •FAVEQf 1000) iFRENOf 1000) tFEPSi 1000) tFRE JO) 1000) tFSl ( 1000) (FH! ( 100 
1  ) *FH2 ( 100 ) »FH3 ( 100 ) t 6AVEQ( 100 ) *GRENO< 100 ) *GEPS ( 100) »6RE JQ( 100 ) * 

1  GS I ( 1 00 ) *  FGl ( 100 ) *  FG2 ( 1 00 ) *  FG3 ( 100 ) *H£PS ( 1 00 ) tHAVEQ ( 1 00 ) » 

IHRENOC 100) tHREJO) 100) •HSK 100) »X( 30*1) 

FREQUENCY  290(0tl«S ) *305(9) * 540(9) tSASC 2* 1*7 ) *565(9) •56S(2*1*7 ) • 

1  595(9) *597 (5* 1*5) *605(5* 1*5) *635 (0*0*1) *655 (300) *657(2*1*7)* 

1  675(300) *677( 2 *1 *7) *696(300) *697 (5*1 *5) *712 (0*0*1) *716(50) *717 ( 

1  2*1*7) *735(50) *737(2*1*7) *755( 25 ) *756( 50 ) *757 ( 5*1 *5 ) *822 (0*1*10) # 
1  a27(0*l*7)*600(300)*792(0*l*3) 

COMMON  Hl*H2*H3*EPS*AVEQ*REJO*RENO 

COMMON  FH1*FH2»X*H1H*FH3*GAVE0*6REN0*GEPS*GREJ0*GS1*FG1 
COMMON  Fl*Sl*F;!*F3*FAVCa*FREN0*FEPS*FRCJ0*FSI 
COMMON  FG3*HEP>*HAVEa*HRENO*HREJO*HlL*FG2*HSI 
NOZ  >  0 

300  read  tape  7*  V1*V2*V3*CE*PE 
READ  TAPE  7  •SL*NR*TL 

READ  TAPE  7  •(  H1(1|*EP  *EP*EP*EP*EP*EP*EP*EP*  I  ■  l*NRt 
HlL(l)  -  HKD 
HlMd)  -  HKMR) 

IF(NS-1)319*315*305 
305  00  310  J  >  2*NS 

READ  TAPE  7  •  V1*V2*V3*C£*PE 
READ  TAPE  7  •  SL*NR*TL 

READ  TAPE  7  •(  HI ( 1 )*EP*EP*EP*EP*EP*EP*EP*EP* 1  >  1*NR) 

HIU(J)  -  Hid) 

310  HIH(J)  ■  HKNR) 

315  REWIND  7 

SMHl  >  Hlt(l) 

560  00  560  1  •  1*NS 

565  lF(Hlt(l)  -  SMHl)  550*560*560 

550  SMHl  ■  Hlt(I) 

560  CONTINUE 

MH1>INTG(SMH1) 

BIGHI  ■  HlH(l) 

565  00  560  I  >  1*NS 

566  IF(BIGH1  •  H1H( I ) )570*580»560 
570  BlGHl  -  HIH(I) 

560  CONTINUE 

IGH1«INTG(BIGH1) 

IJ>1 
590  L>0 

IF(IJ/2)100*100*110 
100  IJ-2 

GO  TO  595 
110  IJ-1 

595  DO  630  (■1*NS 

GO  T0(5«6*596)*1J 

596  READ  TAPE  6  *  V1*V2*V9*CE*  PE 
READ  TAPE  6  *  Sk*NR*TL 

READ  TAPE  6*(  HK 1 1 *H2( 1 ) •H3( I ) vEPS) I ) vAVEOI 1 1 vREJOI 1 1 *REN0( 1 ) * 

1  PL  aKl)*  I  ■  1*NR) 
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GO  TO  597 

598  READ  TAPE  8  *  V1»V2*V3»CE»  PE 
READ  TAPE  8  t  SL*NR*TL 

READ  TAPE  8*(  HI ( 1) tHS ( 1 ) «H3 ( 1 1 tEPSi I ) •AVEQl I ) tRE JOI 1 ) tRENOI 1) • 
I  PL  •SKI)*  I  >  1*NR) 

597  1F(INTG(HK1)  I-MH1)600*600*630 
600  DO  620  J  •  1*NR 
605  IF(  INTGIHKJ)  )>MH1)620*610*630 
610  L  ■  L  •  1 

FKL)  >  HKJ) 

F2(L)  >  H2(J) 

F3(L)  -  H3(J) 

FAVEO(L)  >  AVEQ(J) 

FRENO(L)  ■  RENOIJ) 

FEPS(L)  >  EPS(J) 

FREJOtL)  ■  REJQ(J) 

FSKL)  ■  SKJ) 

620  CONTINUE 
630  CONTINUE 
REWIND  6 
REWIND  8 

635  IF(L)825*825*650 
650  SMH2  «  F2I1) 

655  DO  670  I  altL 

657  IF(F2(II  -  SMH2I660*670«670 

660  SMH2  ■  F2(I) 

670  CONTINUE 

MH2aINT6(SMH2) 

8I6H2  -  F2(l) 

675  DO  690  1  >1(L 

677  1F(BIGH2  -  F2 1 1) )680*690*690 

680  BI6H2  •  F2(n 

690  CONTINUE 

IGH2>INTG(B1GH2) 

695  LI  >  0 

696  DO  710  1  >  1»L 

697  1F( 1NTG(F2(1I )-MH2)710»700*710 
.700  LI  >  LI  ♦  1 

FHKLl)  ■  FKII 
FH2IL1)  •  F2(ll 
FH3IL1)  •  F3(II 
GAVEO(Ll)  a  FAVEQ(l) 

GRENOILl)  a  FRENO(I) 

GEPSILl)  a  FEPSdl 
GREJOiLlI  a  FREJQUI 
GSKLl)  a  FSIdl 
710  CONTINUE 
712  IF(L1)820(020*715 

715  SMH3  a  FH3(1I 

716  00  730  1  a  IfLl 

717  IF(FH3dl  -  SMH3)  720*730*730 
720  SMH3  a  FH3d) 

730  CONTINUE 

MH3alNTG(SHH3) 
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B1GH3  >  FH3(1I 

739  00  790  I  >  liLl 

737  1F(BI6H3  -  FH3<Ilt  740>790»790 

740  B16H3  >  FH3n) 

790  CONTINUE 

lGH3aINTG(B16H3) 

1B1G*0 

MH3  ■  HH3  >  1 
JBlGa  1GH3-  MH3 
799  DO  800  I  ■  ItJBlG 
MH3  >  MH3  4  1 
L2-0 
Z1P>0.0 
ZAP>0*0 
ZOP>0«0 

796  00  790  J  ■  1«L1 

797  1F(MH3-1NTG(FH9(JII I790*760*790 
760  L2-1 

M1P>FH1(JI 

WAPaFHZtJI 

WOP>FH3U> 

WUP-GEPS(J) 

Z 1 P>Z 1 P4GR ENO ( J I •GAVE04J ) 

ZAP>ZAP4GRENO(Ji 

ZOP-ZOP«GREJO(J> 

ZUP-GSKJ) 

790  CONTINUE 

792  IP (L2) 600 *8 00 *799 

799  1G1G-1B1641 
FGKlBlQt-VlP 
FG2I1B16)>WAP 
FG3(1BIG}>W0P 
HAVEQ(1BIGI«Z1P/ZAP 
HREJ0(1B1G}>Z0P 
HSKIBIGI-ZUP 
HEPS(1B1GI«WUP 
HRENO(lBiGI>ZAP 

800  CONTINUE 

810  NOZ  >  NOZ  4  1 

WRITE  TAPE  7f  I  BIG 

WRITE  TAPE  7t(  FGl ( 1) •F62t 1) »FG9I U tHAVEOI I ) •HEPSII I tHREJOl 1 1 • 
1  HRENO(I)(  HSKD*  I  >  ItlBlG) 

820  NH2  ■  MHZ  4  1 

822  IF(  IGH2  -  MHZ  1829 •699t69? 

829  MHl  >  MHl  41 

827  IF<  I6H1  -  MH1)830*990»990 

830  REWIND  7 
REWIND  6 
REWIND  8 
RETURN 

ENOI0fl»0*0*0l 
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SUBROUTINE  A!?.SI  (  SYAM  i  AS  »BS»CS  •  ASl  1 6S1  iNASP  •  AN  *S21  iNOZiEXtSFtNI  • 

1  BTlttlT2> 

DIMENSION  SXT( 100) .SIG( 100) .YSN< 100) »YSD( 100) »WI ( 100) *AVEO( 1000) i 
1  EPS(lC00)i  REJO(  1000)  tSI  ( 1000)  lYST  ( 100)  •AN(20)  lASI  2.0)  * 

1  ASl  l.'ij  )  i.tSI  2C  )  tbSl  120)  tCSI  20)  iHl  (1000  )  iHZI  1000)  tHS  ( 1000)  • 

1  REfJOt 1000) *ZXY(350 ) .X (50»1) »PL1 1000 )• TXT ( 100 ) »AT1(19)»AT2(19) 

1  tSAOdSi) 

COMMON  HliH,2iH3  »EPS  i  AVEOiMEJCltkENO 

COMMON  SXTtYST.XtWI  i YSN » YSO t S IG *2XV »PL tS I  »TXT »AT 1  * AT2 *SAG 
FREQUENCY  10 ( 40 ) lEO ( 4 ) 1 32 ( A  0 ) «40 ( 15 ) *47 ( 300 ) *60 ( 40 ) *  100 ( 40 ) * 

1  1151 15) til? (300) *48(0 il5«l) *50(011*10) (SSI 20il *4) 

SIMAX  »  (l.O/SYAM  )  »»  2 
NI  «XINTF(  SIMAX/.S2I) 

.SXT(  1  )«  S2I  /  2.0 
10  DO  30  I  =  2*NI 

30  SXT(  I  )■  SXT(  I-l)  +  .S2  1 
DO  31  I-1*NI 

SIG(I)  ■  0.0 
20  DO  31  J  ■  1*NASP 

31  SIG(I)  «  SIG(l)  +  AN(J)*(  AS(J)*EXPF(-AS1(J)*SXT(  I  )  )♦ 

1  bS(  J)  *EXPF{-  3S1(  J  )*.*XT  (  I  )  )  ■*-CS(J))**2 

32  DO  35  I  ■1*N1 
YSN(l)  a  0.0 
YSU(I)  ■  0.0 

35  wl ( I )  a  0.0 

40  DO  95  NIP  a  1*N0Z 

45  RFAD  TAPE  7*  IBIG 

READ  TAPE  7  *  ( ABl  .Bi.U  « CB 1  *  A VEO  (  I  )  *  EPS  (  I  )  *  RF  JO  (  1  )  *RENO  (  I  )  * 

1  SKI  )*  la  l*l.:iIG) 

47  DO  90  I  a  l.IHIG 

SIX  a  0.0 

SO  IF(AVEQ( I ) )4e*4e*b0 

48  IF(RF.JQ(  I  )  )60*60*90 
60  DO  80  J  a  1,NI 

SIX  a  SIX  ♦  S2I 
65  IF(SIX  -  SI ( I ) )  80*  70*70 
70  YSN(J)  a  EPS(I)  «  SIG(J)  ♦  YSN(J) 

YSD( J)  a  YSD( J)  ♦  AVEO( I  ) 

WI(  J)  a  wK  J)  4-  1.0 
GO  TO  90 
80  CONTINUE 
90  CONTINUE 

95  CONTINUE 
REWIN(3  7 

00  90  Ial*NI 

IF( YS0( I )-1.0E-l0 196*96 *97 

96  YST( I )a0.0 
GO  TO  98 

97  YST(I)  a  L0GF(  YSN(n/  YSO(I)) 

98  CONTINUF 

EXa.5 

DO  111  Jal*19 
WXX  a  0*0 
WY  a  0*0 


48 


NAVAL  liSIAICH  IA80IAT0IV 


WX  ■  0*0 
WXY  >0t0 
WT  B  0*0 
EXBfX+.5 

100  00  110  I  ■  1*N1 

TXT(  n».SXT(  I  )*»EX 
WT  B  WT  ♦  WI ( I  1 

WXX  B  WXX  +  WKIJ  »  TXT(I)  •  TXT(I) 

WY  B  WY  ♦  wn  n  *  YSTt  I  ) 

WX  B  WX  +  WI ( I  I  *  TXT< I » 

WXY  B  WXY  +  WI(I)  *  TXTin  ♦  YST(I) 
lie  CONTINUE 

AT1(J)b(  WT*  WXY  -  WX  *  WY  )/(WT  *  WXX  -WX**2 ) 

AT  2 ( J 1 B ( WY-AT 1 ( J ) *WX ) /WT 
SAG( J)B0t0 
DO  lllllBltNl 

IFIYSOI  I  )-1.0E-lO»ini.llll»1110 
nio  SAG(  JlB.SAfil  J  )  +  (  AT2«  J)+AT1(  J)*TXT  t  D-YSTI  I  I  )*»2 
nil  CONTI  Ml  IP 
1.11  CONTINUE 

.<;.SAGB.<;Ar,(  i ) 

00  113  IbU19 

IF^(SAQ(  1  )-1.0E-10)113*113*12 
1?  IF(SSAG-SAG(  mil3*112»ll2 

112  tlTlBATl  (  I  ) 

;3T2bAT2  (  I  ) 

EXB0.5+0.5»FL0ATFt 1  ) 

SSAGbSAGI 1 ) 

113  CONTINUE 

SF  B  EXPF(«T2) 

11300120  J  B 1 ,  N02 
READ  TAPE  7  •  IBIG 

READ  TAPE  7  •  (HI  ( 1)  »H2  (  1)  *H3  (  I  )  •  AVr»i  (  I  )  *EPS  (  I  )  tREJOI  1)  tRENOI  I )  * 
1  SI ( 1 ) •  1  B  ItlBlGI 

117  00  130  1  >  ItIBIG 

130  AVEG(  I  )  B  AVEOd  1  «  SF 
WRITE  TAPE  0.  IBIG 

120  WRITE  TAPE  8 t ( HI ( 1 1 tH2 ( I > »H3 ( I ) • AVEO ( 1 1 • EPS ( I ) tREUGI I ) iRENOI 1) • 
1  SI ( I ) t  I  B  1»I3IG) 

REWIND  7 
REWIND  8 
RETURN 

ENO(2«2*2t2t2) 
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SUBROUTINE  0UTPUT(N0Z*NS(N1»BT1*BT2»SF*EX»N01BI 
DIMENSION  HI ( 1 000 i t H2 ( 1 000 ) tHS ( 1 000 ) • AVEQ Cl 000 1 tEPS 1 1000 ) •RE JO 
1C  1000) tRENOClOOOi tPt  C 1000) tSl C 1000) tRAVEOC 1000) •WEI C 1000) •VI ISO) • 
2V2CSO>^Y3ISO)^M1I1000)^M2I1000UMSI1000)^YAVEOC50)^YRAVEO(SO)^ 
SMRENOC 1000 ) • YWE 1 C  50 1 •NRE JO  C 1000) •MEPS 1 1000 ) fSXT  C 100 ) •VST  1 100 ) • 
AXCSO)^  BCBI450)^YS1 ISO)^WIC100) 
COMMONH1^H2^HS^EPS^AVEO^REJQ^RENO^SXT«YST^X^WU  YS^VAVEO^YRAVCO 
l^VWEl^BCB^PL^St  •RAVEQ^WEUM1^N2^MS^MREJ0^NREN0^MEPS^Y1^Y2 
FREQUENCY  190( IS) •19SC2S) •2«0C0^1^10) •SOOl SO^l *0) •940C0^1>20) • 
1362CS)^3B3l900)^S90C0^1^10)^499t0»3^1)^4A0l0^1*9)*S0SCS*l*S0)* 
24B0C2^1^1S) 

ON-1.0 

N-O 

NOIB-0 

NNNbI 

190  00  330  J-l^NOZ 
READ  TAPE  B^IBIG 

READ  TAPE  B •CHI C 1 1 •H2 C 1 1 •H31 1 ) •AVEQI 1 ) tEPSC 1 ) •REJOt 1 1 •RENOC 1 ) • 
1S1CI)(1>1^IBIG) 

N0IBaN0IB4lBlG 
199  00  320  1-l^IBlG 
220  RAVEOC1)>SORTFIAVCO(1)) 

S1C1)>SQRTFCS1C1)) 

260  1FCREN0II))270^270^280 
270  WElCll-0«0 
GO  TO  290 

2B0  WEIClI-RAVEQCD/SQRTFCRENOCIt) 

290  N>N^1 

Y1IN)-H1CI) 

Y2CN)>H2C1 ) 

Y3CN)>H3Ct) 

M1CN)«1NTG<H1CI)) 

M2CN)>1NT0CH2CI)) 

M3CN)>1NT6CH9C1)) 

299  VS1CN)«S1C1) 

YAVEQCN)aAVEOCl) 

YRAVE0CN)>RAVEQCI) 

MRENOCN)>RENOC 1 I^.OOOl 
YWEICN)-WETI1) 

MREJ0CN)>1NTGCREJQC I )) 

MEPSCNr«INT6CEPSCl) I 

300  IFCN-29)320^310^310 

310  IFCSENSESWITCH  6)311«312 

311  PRINT  S40 
PRINT  600 

PRINT  S60^4M1CK)^M2CK)^M3CK)^VS1IK)^YAVEOCK)^YRAVEOCK)^MEPSIK)^ 
1MREN0IK)^MREJ0IK)^K-1^29) 

312  WRITE  OUTPUT  TAPE  9^  S40 
WRITE  OUTPUT  TAPE  9^600 

WRITE  OUTPUT  TAPE  9»  S60»  (MICK) »M2 (H) .M3(K ) »YSI (K) *YAVEO(IC) »YRAVE 
10 I K ) •MEPSC K ) tMRENOC K ) •MREJOC K ) •K«l *29 ) 

WRITE  OUTPUT  TAPE  9 •SBO^ C VICK) •YE CK) •V3I K ) t  YAVEOI R ) tYWEl CK) (ON^ 
1VSIIK)*K«1^N) 
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320  CONTINUE 
330  CONTINUE 
340  1FINI360«360*330 

330  IFtSENSESWITCH  4)331*332 

331  FRINT  340 
PRINT  330 

PRINT  340 • ( Ml ( K ) (Ha I K ) VMS I K ) • YS 1 1 K I • YAVEO 1 K ) • YRAVCOI K ) fMCPS ( K ) • 
IMRENO ( K I tNRE JO ( K ) • Ral tN ) 

332  WRITE  OUTPUT  TAPE  9*340 
WRITE  OUTPUT  TAPE  9*330 

WRITE  OUTPUT  TAPE  9*540*(M1 (K) *H2(K) *M3 (K) *YSI (K) *YAVEO(K)  *YRAVEQ 
1 ( K ) •MEPS ( K ) *MRENO ( K I *MR£ JO U ) * K> 1  * N ) 

WRITE  OUTPUT  TAPE  3*3R0*IV1IK)*Y2(K) *Y3(IC) *  YAVEO(K) *YWE1 (K) *ON* 
lYSI(K)*Kal*N) 

340  END  FILE  3 
REWIND  3 
REWIND  t 

IFISENSESWITCH  4)333*341 

333  PRINT  400 
PRINT  490*N01B 
PRINT  330 

PRINT  390*BT2*Bri*CX*SF 
PRINT  400 
PRINT  410 
PRINT  330 

PRINT  420*(SXT(K)*YSTIK)*Wl{K)*Ral*Nl) 

341  WRITE  OUTPUT  TAPE  9*600 
WRITE  OUTPUT  TAPE  9*A90*N01B 
WRITE  OUTPUT  TAPE  9*330 

WRITE  OUTPUT  TAPE  9*390*BT2*BT1*EX»SF 
WRITE  OUTPUT  TAPE  9*400 
WRITE  OUTPUT  TAPE  9*410 
WRITE  OUTPUT  TAPE  9*330 

WRITE  OUTPUT  TAPE  9*420*ISXT(K)*Y$T(K)*WlCK)*Kal*Nl) 

342  DO  820  K>1*NS 

READ  TAPE  6*U1*U2*U3*CE*PE 
READ  tape  4*SL*NR*TL 

READ  TAPE  4*(H1(  J)*H2I  J)  *H3(  J)  *CPSI  J)  *AVEOI  J)  *REJO(  J)*RENOU  * 
1PL(J)*S1(J)*J-1*NR) 

NUI-INTGIUI) 

MU2>1NT6(U2) 

MU3aINT6(U3) 

MSt>lNT6ISL) 

343  DO  430  1>1*NR 
AVEO(ll>AVEO<l)«SF 
RAVEO(I)>SORTF(AVCOI1)) 

Sll 1)>S0RTF(SI(I)) 

390  IFIRENOI I ) )400*400*410 
400  WE1(I)>0«0 
GO  TO  420 

410  WEI ( I )>RAVEO( I )/40RTF(REN0( I ) ) 

420  Mil I)>INTGIH1(1)| 

M2(I)-INT6(H2(I)) 

M3II)-INTGIH3<I)) 
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MREJO(I)*lNTG(REJO( I n 
MRENOCI )>1NTG(RENC( 1 ) ) 

490  MEPS(l)>INTG(EPS(in 

499  ?p[pEU97»437»436 

496  IFiSENSESWITCH  6110*20 
10  PRINT  695 

20  WRITE  OUTPUT  TAPE  9*635 
GO  TO  440 

497  IPISENSESWITCH  6)30*40 
30  PRINT  690 

40  WRITE  OUTPUT  TAPE  9*690 
GO  T0l440*450*469) *NNN 
440  IPICE >460*460*445 
449  NNN>2 

490  IPISENSESWITCH  6)  50*60 

90  PRINT  640*MU1*HU2*MU9*MSL*SP*NR 
60  WRITE  OUTPUT  TAPE  9*640  *HU1*MU2*MU3*MSL*SP*NR 
GO  TO  470 
460  NNN>3 

465  X(K)>X(K)*SF 
IPISENSESWITCH  6)466*467 

466  PRINT  640*NU1*HU2*MU3*MSL*XIK)*NR 

467  WRITE  OUTPUT  TAPE  9*640*MUl*MU2*MU3*M$L*XlK ) *NR 

470  IPISENSESWITCH  6)471*472 

471  PRINT  990 
PRINT  690 
PRINT  600 

472  WRITE  OUTPUT  TAPE  9*  990 
WRITE  OUTPUT  TAPE  9*690 
WRITE  OUTPUT  TAPE  9*600 

400  lPINRl-22 >490*490*800 

490  NR2>NR-NRl-fl 
IPISENSESWITCH  6)  491*492 

491  PRINT  660*IM1(II*M2II  )*M3I  I)*Slin*Pt(!  )*AVEO(  I )  *RAVE01 1  >  *HEPS 
1II)*MREN0II)*MRCUQ( II*I>NR2*NR) 

492  WRITE  OUTPUT  TAPE  9*660* IMII I ) *M2 1 1 ) *H3 1 1 ) *S 1 1 1 ) *PL 1 1 ) *AVE0I I ) *RAV 
lEQI 1 1 *MEPSI I ) *MRENOI I ) *MREUOI I ) * IaNR2  *NR ) 

60  TO  810 

800  NR2«NR-NR1^1 
NR3>NR2'*’23 

IPISENSESWITCH  6)801*802 

801  PRINT  660* IMl 1 1 ) *M2 1 1 ) *H3 1  I >  *SI 1 1 ) *PL 1 1) *AVEOI I ) »RAVEOI 1 ) *MEPS 
1II)*MREN0II)*MREU0( I l*I-NR2*NR3) 

802  WRITE  OUTPUT  TAPE  9 *660* IMII I ) *M2 1  I ) *H3 1  I ) *S I ( I ) *PLI 1 ) *AVEOI I ) *RAV 
lEOI I ) *MEPS 1 1 ) *MRENO 1 1 ) *MREJO 1 1 ) * I>NR2 *NR3 ) 

NRl-NRl-24 

IPISENSESWITCH  6)803*804 

803  PRINT  670  tMUl •MU2 *MU3 *MSL 
PRINT  650 

PRINT  600 

804  WRITE  OUTPUT  TAPE  9*670*MU1*MU2*MU> *MSL 
WRITE  OUTPUT  TAPE  9*690 

WRITE  OUTPUT  TAPE  9*600 
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BOS  lP(Nftl*24)A90*A90»B00 

BIO  WRITE  OUTPUT  TAPE  7*S80»H1 ( 1 ) *H2 1 1 ) »H3 ( 1 ) *  AVEO ( 1 ) *WE1 ( 1 ) *ON •&! ( 1 ) 
WRITE  OUTPUT  TAPE  7*680»4H1 1  1)  *H2  ( 1)  «H3  (  1 1 1  AVEO(  1 1  *WE1  (  1  >  »  SKI)* 
lla2»NRI 
820  CONTINUE 
END  FILE  7 
REWIND  6 
REWIND  7 
RETURN 

SAO  FORNATIlHltlOBH  H  K  L  SINO/L 

I  F**2  F  EPS  W  R) 

950  FORMAT  I XH  ) 

960  FORMAT! XH0«12H  3I4»FX8.7»F19.3*F16.4.3I9) 

930  FORMAT! IHl I 

980  FORMAT!6F9.2»F9«6*9H  R»12X) 

990  FORMAT433H  LOO  K  >  A  4-  B!S**X)  A  >F10«9»1XH  B  -FlO 

1*9*12H  X  >F9*9*20H  SCALE  FACTOR  >FI0*9) 

600  FORMAT !1H0) 

6X0  FORMAT49XH  !SIN!}/L)**2  LOG  K  W) 

620  FORMAT4X2H  FX0*6*eH  F9«S»FX3*0I 

630  FORMAT !XHX*37H  U  V  W  WEISSENBERG) 

660  FORMAT !XH0*3I3«9H  LAYERI3*36H  SCALE  FACTOR 

X  >FX0«5*26H  NO*  OF  REFLECTIONS  >15) 

639  FORMAT  I XHX*36H  U  V  W  PRECESSION) 

660  F0RMAT!XH0*2H  3I4*FX9*7*FXB*6*FX9.3*FXB*6»3IX0) 

670  FORMAT! XHX«3I3«9H  LAYERI3) 

690  FORMAT  !X20H0  H  K  L  SIN!)/L  LP  FACTOR 

X  FM2  F  EPS  W 

2R  I 

680  FORMAT! 9F9«2*FX8*6*9H  R»X2X) 

690  FORMAT !X HO tAOH  TOTAL  NUMBER  OF  INDEPENDENT  REFLECT  IONS  17) 

ENOI2t2f2*2«2) 
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Moow 

OIMFNSION  HI ( 1000) iH2 ( 1 000 ) *H3 ( 1000 ) lAVEOt 1000) lEPS) 1000 ) tREJO 
1  ( 1000) *  RENO ( 1000 ) tPL ( 1000 ) .SI ( 1000) .RAVEOl 1000) .WEI  (  1000 ) lYl  (  50)  . 
2Y2( 50) ,Y3( 50) .Ml ( 1000) .M2( 1000) »M3 ( 1000 ) . YAVEGK  50) .YRAVEOI 00 ) » 
3MREN0( 1000) .YWEI ( 50 ) .MREJQ( 1000) .MEPS( 1000 ) •  BS 1 ( 2 0 ) tCS ( 20 ) * 
AX(50.1 ) .BCB (450) .YSI ( 50 ) »W1 ( 100) »AN( 20) .AS( 20) »AS1 ( 20) »BS( 20) 
FREQUENCY  190(15) .  195 ( 25 ) *260 ( 0 . 1 » 1 0 ) • 300 ( 5 0 . 1  »0 ) ♦ 34 0 ( 0 1 1 *20 ) * 

1362 ( 5 ) .363 ( 300) .390 ( 0 . 1 . 10 ) .435 ( 0 » 3 . 1 ) .205( 5 ) .  805(5.1.30). 

2480(2. 1.15). 364(5) 

180  REWIND  4 

READ  INPUT  TAPE  4.700 
READ  TAPE  4.NA50 

READ  TAPE  4 . ( AS ( I ) . AS  1 ( I ) .RS ( 1 ) .BS 1 ( I )  .  *5 ( I  ) . AN ( I ) . I « 1 .N ASP ) 

READ  TAPE  4  .NOZ  .NS  . N  1  .RT  1  .RT2  . SF  . EX  .NOl  :•) .  NR T 

REWIND  5 

PFWTMD  6 

REWIND  7 

REWIND  8 

WRITE  OUTPUT  TAPE  5.740.N0I0 
WRITE  OUTPUT  TAPE  7.7Afl  .NRT 
QN=1 .0 
NaO 

185  AE2aO.O 
AE2C=0.0 
AE=n.O 
AECaO.O 
AE21 aO.O 
AE21C=0.0 

190  DO  330  Jal.NOZ 
READ  TAPE  8.IHIG 

READ  tape  8 .  (HI (  I  ) .h2 ( I ) .H3( I ) .AV£0(  I  )  .EPS(  I  » .REJ0( 1 ) .REN0( I ) » 

ISI ( I ) .lal.IRIG) 

195  DO  320  lal.IBlG 

IF( AVE0( I )-1.0E-l 0)22 0.220. 200 
200  SlGaO.O 
205  DO  210  L=1.NASP 

210  SIG*SIG*AN(L )*( AS ( L ) *EXPF ( -ASl ( L ) *S 1(1) )  t-BS ( L) *E XPF ( -PS  1 (L)*SI ( 1 ) ) 
14CS(L) )**2 

AVF0(  I  )aAVFO(  I)  *EXPF(RT1*(  SI  (  n#*EX  )  )  /  (  EPS  (  I  )  *S  I G  ) 

220  RAVFO( I )»SORTF(AVEO( I )  ) 

AVEOd  )aAVEQ(  I  )-1.0 
AF.-AE  +  RAVEQ  (  1  ) 

AECaAFC+l .0 
AF21=AF21+AVFO( I ) 

AE21C»AF21C'M  .0 
AE2aAF2+ARSF(AVFO( I  )  ) 

AE2C»AE2C+1.0 
SKI  )aSORTF{SI  (  I  )  ) 

260  IF( RFNO(  I  )  )270.270.280 
270  WEI ( I ) aO.O 
GO  TO  290 

280  WEI ( I )al ,0/SQRTF(RENO( I ) ) 

290  N-N+l 

Yl(N)=Hl (I ) 
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Y2(N)«H2(I  ) 

Y3{N)«H3( I  ) 

M1(N)»INTG(H1(1 »» 

M2(N)»INT6(H2( I  I) 

M3(N)»INTG(H3( I  >  » 

YSI (N)-S1( I ) 

YAVF0(N)»AV?0( I  1 

yravfo(n)»rave:o(I  i 

MRFNn( N ) aRFNO( I )+»0001 
YWEl (N)iwEI (1 ) 

MREJ0(NI=REJ0( I )+.0001 
MEPS(NI»EPS( I  l  +  .OOOl 
300  IF(N-25)320.310.310 
310  WRITE  OUTPUT  TAPE  9»540 

WRITE  OUTPUT  TAPE  9t6O0 

WRITE  OUTPUT  T APE9 . 560  * ( Ml ( <» *  M2 ( Kl .M3 ( K ) t YS I ( < I .YAVEO(K) »YRAVEQ(< 
1  )  tMEPSIK  1  *MRENO(K»  iMREJOU  I  »K*1  *25  1 
WRITE  OUTPUT  TAPE  5  *680.  ( M 1  (  K  )  .  M2  ( K' )  .M3  (  K  )  .  Y  AVEO  (  K  )  t  K»  1 . 2  5  ) 

M  =  f) 

320  CONTINUE 
330  CONTINUF 
340  !F(N)360*360.350 
350  WRITE  OUTPUT  TAPE  9*540 

WRITE  OUTPUT  TAPE  9*550 

WRITE  OUTPUT  T  APE9 . 560  »  (  Ml  (  K  )  .M2  (  IC>  .M3  (  K  )  *  YS  I  (  lO  *  Y  AVEG  (  <)  *YRAVEQ(IC 

1  )  *MEPS(K  )  .MRENOIK  )  •KREJUI  K  I  *i<*l  *N  ) 

WRITE  OUTPUT  TAPE  5  *680  .  ( M 1  (  K  I  .M2  I  I  .M3  t  <)  .  YAVEQ  (  K  )  *IC«  1  *  N  ) 

360  END  FIl.F  5 
RFWINO  5 
RFWIMO  fl 

WRITF  OUTPUT  TAPE  9*530 

WRITF  OIJTPUT  TAPF  9 , 590  ♦  RT2  .  «T  1  *  FX  .  SF 

361  ZAFsAF/AEC 
?AF?sAF2/AF2C 
ZAF21  =AF2l  /APZK 

WRITE  OUTPUT  tape  9*600 

WRITE  OUTPUT  TAPE  9 . 690 *2  At  *  2 AE 2  *  2 AE 2 1 

362  00  820  <=1 *NS 

READ  TAPE  6 *U 1  * U2 .U3 * CE »PE 
READ  TAPE  6*SI.»NR*TL 

READ  TAPE  6 . ( H 1 ( J ) . H2 ( J » »H3 I J ) *F PS  I J ) * AVEO ( J ) *  RE JO  I J ) *  RE  NO ( J ) » 

10L( J)  .SI ( J) .J=1»NR) 

MU1=IMTG<U1 ) 

MU2*IMTri(U2  ) 

MU3=INTG(U3 I 
MSL=INTG(SL ) 

363  00  430  1  =  1  .NR 
SIG=0.0 

370  no  380  L  =  1  .NASP 

380  SIG=S IG+AN ( L )*( AS( L ) •EXPFI -ASl (L)»SI(II)4BS(L)*EXPF(-BS1(L)*SI(I)) 
1+CS(L) )**2 

PL( I )»FXPF(BT1*(SI (I )**EX) ) /(EPS( I '*SIG) 

IF( AVFOI I )-1.0E-10)365.385.381 

381  AVFOI  I  )  =AVFO(  I  »*Pl.<  I  )*.SF 


NAVAL  RiSIARCH  LABORATORY 


55 


I 


RAVPO(  I  )=SORTF(  AVFn(  1  M 
I  AVF0( ! ) »AVFO( I  1-1 

!  SI ( ! )»FORTF(St ( I ) ) 

»  390  IF(RFN0( I ) 1400.400*410 

'  400  WFI ( I  1*0.0 

00  TO  4?0 

410  WEKI  l*1.0/S0RTF(REN0t  I  1  1 
420  Mil  1  1»IMTG(H1(I  1  ) 

M2( I 1=INTG(H2 (Ill 
M3( t )=INTG(H3( T 1 1 

429  MRFJOI I l»RFjn( I l+.OOOl 
WRENOt I l*RFNO( T l+.OOOl 

430  MEPS(  I  1  *EP.S  (  I  l  +  .OOOl 
NR IsNR 

433  IFIDF )437«437.436 

436  WRITE  OUTPUT  TAPE  9.635 
GO  TO  450 

437  WRITE  OUTPUT  TAPE  9.630 

450  WRITE  OUTPUT  TAPE  9 . 640 .MUl »MU2 ♦MU3 .MSL »NR 

470  WRITE  OUTPUT  TAPE  9.550 

WRITE  OUTPUT  TAPE  9.650 

WRITE  OUTPUT  TAPE  9.600 

480  IF( NRl-22 1 A90.490.fl00 
490  NR2=NP-NP1+] 

WRITE  OUTPUT  TAPE  9 .660 . ( Ml ( ! 1 .M2  11 1 .M3  I  I  1 . S I ( I  1 ♦ PL ( I  1 . AVEO 11 1 » 
IRAVFOI  I  1  .MFPSI  I  1  .MRFNOI I  I  .“REJOI  I  1  .I.»MR2.NP  1 
GO  TO  810 
800  NR2*MR-NRl4’l 
NR3=NR2+23 

WRITE  OUTPUT  TAPE  9 . 660 . ( MI ( I  1 .M2 t 1  1 .M3  I  I  1 .S I ( II . PL ( I  1 . AVEO I  1 1 . 
IRAVEOI I  1 .MEPS( I  1 .MRENOI 1 1 .MREJOI I  1 . I*NR2*NR3 1 
NRlaNRl-24 

WRITE  OUTPUT  TAPE  9 .670 .MUl *MU2 .MU3 .MSL 
WRITE  OUTPUT  TAPE  9.650 
WRITE  OUTPUT  TAPF  9.600 
805  IF(nr1-?41490.490.800 

810  WRITE  OUTPUT  TAPE  7 .680. I  MU  I  1 .M2 ( I  1 .M3 ( I  1 . AVEO ( I  1 . I *1 .NR  1 
820  CONTINUE 

END  file  7 
REWIND  6 
rewind  7 
PRINT  730 
PAUSE 

540  FORMAT! IHl .106H  H  K  L  SINIl/L 

1  E**2  -  1  E  EPS  W  R1 

550  FORMATIIH  1 

560  FORMAT! 1H0.12H  3 1 4 .F 1 8 . 7 . F 1 9. 5 . F 1 fl , 5 . 3  I  9 1 

530  FORMAT! IHI  ) 

590  F0RmaT!33H  |  OG  5  =  A  ♦  R!S**X!  A  ■F10,5.11H  B  aFlO 

1.5.12H  X  =F9.5.2nH  SCALE  FACTOR  =F10.5) 

600  FORMAT!1H01 

630  FORMAT! 1H1.37H  U  V  W  WEISSENBERGl 

640  FORMAT!1H0.3I3.9H  LAYERI3.36H  NO.  OF  REFLECTIONS 

1  *151 
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635  F0RMAT(1H1.36H  U  V  W  PRECESSION) 

660  FORMAT(lHOt2H  3 14* F19.7*F1B.6 *F 19. 5 tF18* 5 f 3  1 10 ) 


670  F0RMAT(1H1.3I3»9H  LAVERI3) 

650  FORMAT  (120H0  H  L  SINO/L  K 

1  E**2  -  1  E  EPS  W 

2R  ) 

680  FORMAT(3I4.F20.5) 

690  FORHAT(1HO*12H  AVE.  E  >  F10.6*24H  AVE.  /E**2  -  1/  ■F10.6« 

122H  AVF.  E**2  -  1  -FIO.S) 

700  FORMAT  I 72H 

1  ) 

710  FORMAT! I4tF5.4l 
720  F0RMAT(5F8.4.F4.0) 

730  F0RMAT(80H  JOB  FINISHED.  REMOVE  AND  SAVE  TAPE  7.  USE  TAPE  5  FOR 
1  INPUT  TO  JOB  r  ORDERING.  ///////////  ) 

740  FORMAT! 17) 

EN0!0«lt0.0.1) 
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E  ORDER 

DIMENSION  Ml  12400) I  M2 (2400) (MS (2400) (AVI  2400) »K 1(2400) fK2(2400) tK 
1(2400) tQAV( 2400) iR) 2400) 

COMMON  MliM2*M3*AViR 

FREQUENCY  5 ( 5* 1  iS ) 1 30 ( 1200 ) tSO) S *  1  * S ) iB0( 600 ) •901 S 1 1 *5 ) * 
1170(lil«3)tie0(25) 

REWIND  5 

READ  INPUT  TAPE  SilOOltNOIB 

1  IF(SOOO*N01B)2*3i3 

2  NIBO-SOOO 
NOIB-NOIB-3000 
NBI0>1 

GO  TO  4 

3  NBIO-2 
NIB0>N01B 

4  K-O 
L«0 

DO  20  I«1*N1B0 

READ  INPUT  TAPE  S* lOOOtNl •N2«N3*0 

5  IF(0)20«10«10 
10  L«(.+  l 

Ml(L)*Nl 
M2(L)«N2 
M3(L)»N3 
AV(L)>Q 
20  CONTINUE 
30  DO  160  J>1«L 
40  DO  70  K«1*L 
30  1F(AV(K))70*60«60 
60  BAV«AV(K) 

GO  TO  60 
70  CONTINUE 
SO  DO  110  1>K*L 
90  1F(AV( I)-aAV)ilOflOO«100 
100  BAV>AV(I) 

11«I 

110  CONTINUE 

Kl( J)«M1(I1  ) 

K2( J)-M2(I1  ) 

K3( J)»M3( 1 1  ) 

OAV(J)aAV(  ID 
R(J)«SQRTF(AV( Il)+1«0) 

AV(I1)>-1«0 
160  CONTINUE 

WRITE  OUTPUT  TAPE  9*1006 
WRITE  OUTPUT  TAPE  9*1001*NOIB 
MaL/50 

170  IF(M-l)200tl80.180 
160  DO  190  lal*M 
N-l»50-49 
Nl«N+24 

WRITE  OUTPUT  TAPE  9.1002 
WRITE  OUTPUT  TAPE  9*1003 

190  WRITE  OUTPUT  TAPE  9 *1004  * ( K1 ( J) *K2 ( J ) *K3 ( J ) *OAV ( J ) *K ( J ) *K1 ( J425 ) * 


i 
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1K2(  )  iK3  ( J-f25  )  *OAV(  •!%  ( J-f2S  )  I  J>NtNl ) 

200  WRITE  OUTPUT  TAPE  9tl002 
WRITE  OUTPUT  TAPE  9*1003 
L1»L-50»M 
L2-L1-25 

IF(L2)220i220*210 
210  N»M»50+1 
N1-N+L2-1 

WRITE  OUTPUT  TAPE  9 *1004*  (  Kl  (  J )  *<2  (  J )  (KS  (  J)  *0AV(  J )  *R  ( J )  *K1  (  J<f2S  )  * 
1K2(  J-f23  )  *<3  ( J-^25  >  *QAV(  J-f25 )  tR  ( J4-23  )  •  J>N*N1 ) 

N2-N1+1 

N3>N-*-24 

WRITE  OUTPUT  TAPE  9 • 1003 » ( K1 ( J) *<2 ( J ) •K3 ( J) *0AV( J) *R ( J ) • J-N2 •N3 ) 
GO  TO  230 
220  N»M*50+1 
Nl»M*50+Ll 

WRITE  OUTPUT  TAPE  9*1003  •(  K1  (  J )  tlL2  ( J )  *1(3  (  J)  *QAVi  J)  tR  (J  )  tJ>N*Nl ) 
230  call  OUI (K1 *K2*K3*0AV*L> 

GO  TO( 1*2401 *NBIO 
240  PRINT  1007 
STOP 

1000  F0RMAT(3I4*F20*St 

1001  FORMAT! 17) 

1002  FORMAT! 1H1*107h  H  K  L  E»*2  -  1  E 

1  M  K  L  E«*2  -  1  E) 

1003  FORMAT!1HO) 

1004  FORMAT!lH0*113*2i4*2F13*9*122*214*2F13*3) 

1005  FORMAT!lHO*113*2i4*2F13*3) 

1006  FORMAT !1H1) 

1007  FORMAT  1 1H1*28H  JOa  FINISHED*  SAVE  TAPE  3*/////) 

EN0!1*1*0*0*1) 
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SUBROUTINE  OUl IKl iK2tK3*QAV*LI 

DIMENSION  K1 1 2400 1 tK2 ( 2400 1 •KSI 2400 ) tOAVI 2400 ) f KU I S00*i ) (LUI 300* 
lS)*MUI300tO)*Q(300*OI  •S(300»0I 
COMMON  RUtLUfHUtQ  tS 
J1>0 
J2>0 
JS«0 
J4«0 
JS-O 
JfO 
JTaO 
JS«0 

10  DO  300  J>1(L 

20  IF(Kl(JIH(lUt/2H(l(J)/2|40*30*40 
30  tl*2 
00  TO  30 
40  I 1-1 

SO  lF(K2(J)-K2(JI/2-K2(JI/2)T0»6O*70 
SO  12-2 
GO  TO  60 
70  12-1 

80  lF(K3IJi-K3(Jt/2-K3(JI/2l95*«0*«S 
90  13-2 

GO  TO  99 
99  13-1 

99  GO  T0(100*200)*I1 

100  GO  T0<110*120)*12 

110  60  TOail*112i«13 

120  GO  T0(121*122)«I3 

200  60  T0(210*220)*12 

210  GO  T0(211*212I*I3 

220  GO  T0(221*222I«I3 

111  J0-J8«l 
RUIJ8«6)-K1(JI 
LU(J8t8)-K2(J) 

MU(J6*6)-K3<J) 

0(J6*0)-QAV(JI 

SU8*6l-SaRTF(QAV(J)4l.0) 

60  TO  300 

112  J9-J9-»l 
IIU(J9*S)-K1(J) 

LUU9t9l-R2(J) 

MU(J9t9l>K3(J) 

Q(J5*3)-OAV(J) 

S<J9»SI>SaRTF(QAV(J)-M«OI 
GO  TO  300 

121  JS-JS-M 
0(JS»6I-0AV(JI 
ILUIJ6*6)-K1(J) 

LU(J6fS)-K2(J) 

MU(J6tS)-K3(JI 

S(J6*6)-S0RTF(QAV(J|41«0) 

GO  TO  300 

122  J2-J241 
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KU( J2*2)>K1(JI 
LU(J2*2)-K2(JI 
MU(J2*2)-K3(J) 

Q(J2*2>>QAV(J) 

S(J2t2)>SQRTF(0AV(JKl«0) 

60  TO  300 

211  J7"J7*1 
KUIJ7*7)-K1(J) 

LU(J7t7)-K2(J) 

MU(J7t7l-K3(J) 

Q(J7*7)-OAVUI 

S(J7*7)-S0RTP(0AV(J)-fl*0) 

GO  TO  300 

212  J3>J3*1 
RU<J3*3)-K1(J) 

LU( J3t3|aK2( J) 

MU( J3*3)-K3(J) 

0(J3«3)-0AV(J) 

$(J3*3>-SORTF(OAV(J)^1*0) 

60  TO  300 

221  J4-J44-1 
KU<J4*4)>K1(J) 

LU(J4«4)>K2(J) 

MU(J4*4)-K3(J) 

Q(J4*4laOAV( J) 

S(J4*4)>SQRTF(OAV(J)^1*0) 

60  TO  300 

222  Jl-Jl-t-l 
RU(Jl*l)<iRl(J) 

LM(J1«1I«K2(J} 

MUU1*1}-K3(J» 

Q(Jl<*ll«OAVUi 

SUl*l)>SQRTF(OAVCJ)4l*0) 

300  C0NTINU6 
LL-Jl 

00  500  Kal>a 
M>LL/50 

310  IF(I4-1I340»320»320 

320  00  330  I-ltM 
N« I *50-49 
N1>N*24 

WRITE  OUTPUT  TAPE  9*1002 

WRITE  OUTPUT  TAPE  9*1003 

330  WRITE  OUTPUT  TAPE  9*1004*(KU(  J*K>  *LU(  J*K)  *MUU*K)  *OIJ*R)  *61  J*K|  • 

1KU( J*2S*K)*LU(J*25*K|*MUI J*29*R) »0IJ*25*K)*S{ J*25»Ki •U*N*N1) 

340  WRITE  OUTPUT  TAPE  9*1002 

WRITE  OUTPUT  TAPE  9*1003 

l.l«LL-90*W 
L2-L1-25 

350  IF(L2)370*370.360 

340  N-M*50*l 
N1»N*L2-1 

WRITE  OUTPUT  TAPE  9*1004*tKUIJ*K) *L.U( J*K> *MU( J*R) *Q( J*K) *SI J*R) * 
1KU( J*25*K)*LU( J*29*K)*MUIJ<*25*KI*OIJ*25*K)*S(J*29*K) *J*N*N1) 
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N2-N14-1 

N3-N<f24 

WRITE  OUTPUT  TAPE  «»100S«(KU( J*K ) •LU( J*K| »MU( J*K) (Ot J*K) tSI J»K) t 
1J>N2»N3) 

60  TO  380 
370  N«M»50+1 
Nl«M»50*Ll 

WRITE  OUTPUT  TAPE  9 1 1003*  I  KUt  J*R}  •LUt  J*K »  *MUI  J»K  t  •0IJ*K)  tBUtK)  * 
IJaNtNl) 

380  GO  TO(390*400t410f420*430(440*4»0»900|»K 
390  LL>J2 

GO  TO  300 
400  LL-J3 

60  TO  300 
410  LL>J4 

GO  TO  300 
420  LL-J3 

60  TO  300 
430  LL-J6 

GO  TO  300 
440  LL>J7 

GO  TO  300 
430  LL>J8 
300  CONTINUE 
RETURN 

1002  FORMAT! IHltlOTH 

1  H 

1003  FORMAT! IHO) 

1004  FORMAT!1HO«I13*2I4*2F13«3*122*2I4»2F13*3> 

1003  FORMATI1HO«I13*2I4*2F13*3I 

EN0!2«2*2*2*2) 


H  K  L  E**2  -  1 

K  L  E»*2  -  1  El 
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FUNCTION  INTG(A) 

IF(A)20»10t20 
10  INTG-0 
GO  TO  30 

20  INTG>XSIGNF(XINTF(AB5F(AK*1I  »XINTF(AI  ) 
30  RETURN 

ENOU.l,  2 .7 
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